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Abstract: 


Generation  of  induced  seismicity  in  underground  mines  is  controlled  by  various  local 
and  mine-wide  factors,  and  the  seismic  series  are  usually  composed  of  events  of  different 
origins.  Although  it  is  well  known  that  seismicity  in  mines  is  a  multimodal  phenomenon,  it 
is  not  easy  to  identify  its  various  components.  A  set  of  statistical  methods  to  study 
complexity  of  seismic  series  is  proposed.  These  are  the  descriptive  statistics  analysis  to 
recognize  basic  features  of  the  seismic  series,  the  extreme  value  method  to  identify  modes 
of  the  energy  distribution  of  events,  to  estimat  the  limiting  energy  values  for  these  modes 
and  to  distinguish  those  event  generating  processes  which  have  well  separated  modes  of 
energy  and  the  analysis  of  deflection  to  identify  generating  processes  that  are  responsible  for 
certain  dominant  directions  in  the  seismic  series.  The  dominant  directions  of  the 
distribution  of  epicenters  are  described  by  the  deflections  of  straight  lines  coimecting  the 
epicenters  of  every  two  consecutive  events  in  a  series,  measured  from  the  NS  direction.  The 
deflection  then  forms  another  set  of  time  series.  The  series  of  deflections  are  analyzed  by  a 
nonparametric  kernel  estimation  of  the  deflection  probability  density  function;  a  novel, 
objective  quantitative  approach. 

The  proposed  methods  were  applied  to  the  sets  of  events  from  Wujek  coal  mine,  that 
occurred  from  1988  to  1993  in  different  areas  of  the  mine.  Altogether  14297  events  was 
analyzed.  Their  distribution  is  not  random  both  in  space  and  in  energy.  The  analyses 
resulted  in  identification  of  likely  processes  of  seismic  event  generation.  The  analysis  of 
deflection  applied  in  the  pseudo-multivariate  way  can  be  used  to  separate,  from  the  seismic 
series,  events  originated  by  these  identified  processes. 
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1. 


Introduction 


Generation  of  induced  seismicity  in  underground  mines  is  controlled  by  various  local  and 
mine-wide  factors.  Recorded  series  contain  events  originated  by  different  processes  of  rockmass 
fracturing.  Recognition  of  the  complex  structure  of  the  series  and  separating  of  their  components 
can  be  helpful  for  determining  mechanisms  governing  the  seismic  event  generation. 

The  research,  made  within  the  framework  of  the  Special  Contract  SPC-95-4029,  aimed 
to  answer  to  the  principal  question:  Is  it  possible,  by  the  statistical  analysis  of  standard  seismic 
catalogs,  to  evidence  the  fact  that  the  seismic  series  contain  complex  mixture  of  outcomes  of 
different  generating  processes?  When  such  the  question  is  positively  answered  then  it  will  be  also 
possible  to  identify  those  processes  and,  after  complementing  the  results  of  the  statistical  analysis 
with  information  from  other  sources,  to  conclude  about  possible  factors  controlling  the  identified 
processes. 

Statistical  inference  applied  to  samples  of  events  described  by  a  number  of  parameters 
means  that  these  parameters  are  regarded  as  random  variables,  having  specific  probability 
distributions.  Their  measured  values  are  supposed  to  be  the  statistical  samples  representing 
features  of  the  populations.  If  the  problem  set  by  us  can  be  resolved  by  statistical  analysis 
methods  then  the  parameters  of  the  seismic  events,  generated  by  difierent  processes,  must  differ 
statistically,  which  means  that  they  must  have  different  probability  distributions.  In  such  the  case 
the  general  probability  distributions  of  parameters  for  the  whole  seismic  series,  which  is  a 
mixture  of  outcomes  of  different  processes,  are  superpositions  of  the  individual  distributions 
related  to  the  particular  processes  and  are,  in  general,  multimodal.  In  a  simplified  but  feasible 
approach  one  can  regard  the  modes  of  these  generalized  distributions  as  identifying  their 
unimodal  components,  hence  as  identifying  different  generating  processes. 

The  above  presented  concept  of  identification  of  different  processes  generating  the 
seismic  events  can  be  realized  in  two  ways.  The  first  one  is  a  univariate  approach  in  which  each 
of  the  parameters  is  studied  individually,  which  leads  to  individual  separating  of  the  modes  of  its 
probability  distribution.  This  approach  is  relatively  simple,  has  a  good  theoretical  backgroimd 
and,  thanks  to  this,  in  many  cases  can  be  used  in  practice.  Its  main  drawback  stems  from  the  fact 
that  no  crossing  relations  between  parameters  are  taken  into  account.  These  crossing  effects  are 
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obviously  present  for  the  obvious  reason  that  the  parameters  describe  the  same  phenomenon. 
The  second  approach  is  multivariate  in  which  the  multidimensional  distribution  of  all  parameters 
is  studied.  It  is  conceptually  more  proper  than  the  first  one,  however  difficulties,  that  appear 
when  attempting  to  determine  exact  characteristics  of  the  multidimensional  distribution  from  the 
sample,  cause  that  this  approach  can  be  used  only  in  a  simplified  way  and  only  when  the  data 
satisfy  very  specific  conditions.  An  important  representation  of  the  univariate  approach  is  the 
cluster  analysis.  This  analysis  has  been  developed,  however,  for  phenomena  described  with 
numerous  random  variables,  which  values  can  be  measured  repeatedly  (Anderberg  1973).  A 
version  of  the  cluster  analysis,  the  single-link  clustering  method  was  used  to  separate 
components  of  the  seismic  series  (e.g.  Frohlich  and  Davis  1986,  Davis  and  Frohlich  1991).  Our 
trials  made  in  this  work  on  the  mining-induced  seismicity  for  the  time  being  disqualified  this 
method  for  our  purposes  when  so  many  events  described  by  so  few  parameters  were  to  be 
studied. 

Within  this  study  we  decided  to  work  out  a  compromise  between  the  univariate  and 
multivariate  approaches,  which  we  have  called  as  the  pseudo-multivariate  approach.  It  relies  on 
repetitions  of  the  univariate  analysis  (the  one  parameter  analysis)  for  subsets  extracted  from  the 
whole  set  of  data  so  that  values  of  other  parameter  (or  parameters)  in  the  subsets  are  kept  in 
selected  intervals.  Comparison  of  features  of  the  analyzed  random  variable  determined  for 
different  subsets  reveals  partially  the  crossing  relations.  Out  of  the  parameters  of  seismic  events, 
which  are  in  our  standard  catalogs  from  mines:  time  of  the  event  occurrence,  coordinates  of  its 
epicenter  and  its  seismic  energy,  the  coordinates  of  epicenter  and  the  energy  have  been  used  to 
split  the  catalog  into  the  subsets.  The  reasons  for  choosing  the  first  discriminator  are  obvious. 
We  expect  that  the  seismic  events  grouped  in  an  area  are  more  structurally  homogeneous  than 
the  events  that  come  from  different  clusters.  The  energy,  from  the  other  hand,  according  to 
world-wide  studies  of  induced  seismicity,  most  probably  discriminates  different  kinds  of  the 
generating  processes.  Two  classes  of  seismic  induced  events  were  distinguished,  in  general  the 
small  events,  that  are  directly  connected  to  mining  works  and,  in  general,  the  strong  ones  which 
are  triggered  by  mining  works  but  which  are  also  controlled  by  tectonic  and  mining  structures 
(Gibowicz  1992,  Johnston  and  Einstein  1992). 
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Statistical  analysis  works  on  statistical  samples  that  is  the  samples  that  are,  simply 
speaking,  constant  in  time.  A  considerable  effort  have  been  done  in  studying  basic  features  of 
our  data  in  order  to  test  their  time  homogeneity.  The  results  are  presented  in  Section  4  of  this 
report.  The  descriptive  statistics  studies  allowed  also  to  determine  global  variability  of  the  event 
parameters,  to  verify  correcmess  of  the  way  of  division  of  the  catalog  into  subsets  and  to 
determine  the  intervals  of  coordinates  and  seismic  energy  for  such  the  division. 

Sometimes  the  problem  of  data  decomposition,  difficult  for  solving  by  methods  directly 
dealing  with  measured  values,  can  be  significantly  simplified  by  means  of  an  indirect  method 
that  uses  such  data  transformation  which  strongly  reduces  the  original  data  set.  This  kind  of 
method  was  in  our  study  the  method  of  the  analysis  of  asymptotic  distributions  of  the  extreme 
values.  Instead  of  studying  of  distributions  of  the  seismic  event  energy  the  distributions  of 
energy  maxima,  evaluated  in  constant  time  intervals,  are  analyzed.  The  method,  presented  in 
Section  5,  was  already  used  in  studies  of  the  seismic  catalogs  from  mines  (Kijko  et  al  1987). 
Here,  however,  is  has  been  applied  for  the  first  time  so  widely  with  the  use  of  statistical  tests  for 
assessing  differences  between  the  empirical  distributions  of  the  energy  maxima. 

On  the  basis  of  numerous  works  describing  both  natural  as  well  as  induced  seismicity  one 
can  expect  that  the  spatial  distribution  of  events  generated  by  a  single  process  will  have  a 
distinct,  more  or  less  linear  direction.  For  the  purpose  of  studies  of  the  directional  character  of 
the  seismic  series  we  have  introduced  a  new  parameter,  called  as  deflection,  constructed  on 
mutual,  bi-directional  migrations  of  two  consecutive  events.  Its  definition  and  its  basic  statistical 
properties  are  given  in  Section  6. 

We  anticipated  initially  that  it  would  be  possible  to  identify  modes  of  the  distribution  of 
deflections  from  the  histograms.  However  subjectivity  of  histograms  causes  that  such  the 
identification  would  have  been  non-unique  and,  of  weakly  reliable.  As  the  things  are  we  decided 
to  find  objective  methods  of  estimation  of  the  probability  density  function  (p.d.f.)  from  the 
sample  data.  The  nonparametric  kernel  estimation  of  p.d.f.  turned  out  to  suit  to  our  purposes. 
This  method,  so  far  as  we  know  for  the  first  time  used  in  the  seismic  series  studies,  is 
insufficiently  developed  from  the  point  of  view  of  specific  practical  applications.  Therefore  we 
had  to  make  much  effort,  to  carry  out  a  lot  of  simulation  and  comparative  experiments  to  be  able 


to  use  it  for  studying  the  multimodality  of  the  distributions  of  deflection.  The  results  of  these 
studies  are  given  in  Section  7. 

In  Section  8  we  present  the  detailed  analysis  of  deflections  for  a  number  of  subsets  of  our 
seismic  catalog,  the  identification  of  generating  processes  of  seismic  events  on  the  basis  of 
modes  of  the  distributions  of  deflection  and  some  attempts  of  relating  the  dominant  trends  of  the 
seismic  series  to  tectonic  and  mining  structures  of  the  mine.  The  method  of  the  directional 
character  of  seismic  series  analysis  with  the  use  of  the  nonparametric  kernel  estimation  of  p.d.f. 
will  be  presented  on  the  4th  International  Symposium  on  Rockburst  and  Seismicity  in  Mines  in 
1997  (Lasocki  et  al.  1996). 

Section  9  contains  an  example  of  separating  of  seismic  events  that  are  outcomes  of  the 
generating  process  identified  by  the  analysis  of  deflections.  Finally,  Section  10  presents 
conclusions  which  we  could  draw  from  the  research  made  within  this  project. 
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2. 


The  area  of  investigatioo 


When  selecting  the  area  for  investigations  we  considered  a  number  of  conditions  that  the 
data  should  satisfy. 

1.  Size  of  sample  sets.  Generation  of  induced  seismicity  in  underground  mines  is 
controlled  by  various  local,  mine- wide  and  regional  factors.  However  the  majority  of  seismic 
events  is  directly  linked  to  the  mining  operation.  The  seismicity  connected  with  mining  can  hide 
effects  of  processes  with  lower  event  rate.  We  wanted  to  test  whether  it  was  possible  to  ascertain 
a  complex  character  of  seismic  catalogs  and  separate  its  components  from  seismic  series  analysis 
and  to  develop  suitable  methods  for  such  separation.  Thus  we  had  to  possess  initially  huge  data 
sets.  In  this  way  the  statistical  inference  could  deliver  certain  conclusions  for  samples  reduced  in 
course  of  separating. 

2.  Completeness  and  uniform  quality  of  seismic  catalogs.  The  requirement  to  have  huge 
data  set  implies  that  the  data  will  come  from  a  considerable  time  period.  For  obvious  reasons 
our  studies  could  be  performed  only  on  complete  and  uniform  quality  catalogs.  Therefore  we 
had  to  be  sure  that,  during  a  period  selected  for  our  study,  all  a  seismological  network  layout, 
quality  of  recording  and  event  parameter  estimation  techniques  remained  unchanged.  These 
were  only  initial  requirements  and  we  also  studied  an  internal  structure  of  our  data  sets  from  a 
point  of  view  of  their  completeness. 

3.  Wide  energy  range  of  the  data.  Polish  standard  seismic  catalogs  in  mines  contain 
limited  number  of  parameters  of  recorded  events.  Each  event  is  described  by  the  time  of  event 
occurrence,  epicenter  coordinates,  some,  usually  poor,  information  about  the  depth  of  focus  and 
the  seismic  energy  at  the  source.  We  expect  that  different  generating  processes  result  in 
differences  in  distributions  of  parameters  of  events.  Following  world-wide  studies  (e.g., 
Gibowicz  1990,  Johnston  and  Einstein  1990),  we  can  suspect  a  particular  significance  of  energy 
distributions.  Triggering  processes  are  supposed  to  generate  events  with  higher  energy  modes 
than  processes  directly  controlled  by  mining.  Therefore  it  was  necessary  to  have  data  covering 
the  widest  possible  energy  range. 

4.  Complex  mining  and  geological  conditions  of  an  area.  We  had  to  have  data  from  an 
area  where  more  than  one  generating  processes  were  likely  to  happen.  Studies  on  structure  of 


induced  seismicity  postulat  that  there  are  at  least  two  different  types  of  this  phenomenon:  the  one 
directiy  connected  with  mining  and  the  other  triggered  by  faults  and  other  weak  zones  (Gibowicz 
1990,  Gibowicz  and  Kijko  1994).  Thus,  to  make  our  study  feasibly,  we  wished  to  have,  in  the 
area  under  study,  structures  that  could  give  rise  to  the  second  type  of  induced  processes.  This 
condition  does  not  mean  that  developed  methods  for  generating  processes  identification  are 
limited  to  the  processes  originated  by  faults  or  similar  structures.  However  we  had  to  have 
chance  to  find  more  than  one  component  in  our  data. 

5.  Availability  of  seismological  and  complementary  data.  In  order  to  form  reasonable 
basis  for  sample  reduction  and  to  validate  conclusions  concerning  existence  of  different 
processes  not  only  seismological  data  were  required  but  we  had  to  have  mining  and  geological 
information  about  the  area. 

The  given  conditions  were  well  satisfied  by  data  fi'om  the  Wujek  coal  mine.  Wujek  coal 
mine  exploits  seams  in  the  central  part  of  the  Upper  Silesian  Coal  Basin  (USCB)  located  in 
southern  Poland.  Structurally  the  area  of  mine  belongs  to  the  main  anticline  of  USCB.  The  area 
is  bounded  firom  the  south  by  double  Klodnicki  fault,  with  some  125  and  120  meter  throw, 
respectively,  from  the  east  by  Wojciech  fault,  with  80  to  120  meter  throw,  and  from  the  west  by 
Srodkowy  fault  with  50  meter  throw.  The  area  is  additionally  fragmented  by  a  number  of 
secondaiy  faults  with  general  the  N-S  or  NE-SW  strike  directions.  According  to  some  recent 
hypotheses  the  rockmass  in  this  part  of  USCB  is  still  not  in  equilibrium  and  the  strongest  mining 
tremors  may  be  resulted  by  a  superposition  of  a  local  stress  field  due  to  exploitation  and  non- 
relaxed  remnant  tectonic  stresses  (Teper  et  al.  1992).  The  faults  present  in  the  mine  area  divide 
it  into  the  following  five,  more  or  less  homogeneous,  parts: 

-  The  north-western  (NW),  the  left  uppermost  part  of  the  main,  roughly  below  the  line 
assigned  by  (17500,  -8500)  and  (19000,  -6000)  coordinates; 

-  The  south-western  (SW),  below  NW  part,  enclosed  from  the  east  by  the  line  assigned 
by  (19500,  -7400)  and  (18100,  -7400)  coordinates; 

-  The  central  (C),  left  to  NW  and  SW  parts,  enclosed  from  the  east  roughly  by  the  line 
of  (19650,  -8850)  and  (17500,  -8700)  coordinates.  Mining  has  been  carried  on  only  in  the 
southern  part  of  this  area; 


10  — 


-  The  south-eastern  (SE),  left  to  C  part,  enclosed  from  the  north  by  the  line  of  some 
(19000,  -10500)  and  (18450,  -8750)  coordinates; 

-  The  north-eastern  (NE),  bounded  by  the  northern  border  of  SE  part  and  the  eastern 
border  of  C  part. 

Because  the  mine  is  partially  just  under  the  town  Katowice  mining  has  been  done  by 
means  of  longwall  technique  in  tabular  excavations  in  some  areas  with  hydraulic  backfill  and  in 
other  by  caving. 

The  mine  is  old  and  there  is  a  lot  of  interacting  remnants  of  past  works  since  1920 
(edges,  pillars  etc.)  at  various  levels.  The  remnants  definitely  influence  stability  of  the  rockmass. 

Both  the  recording  system  and  the  data  collection  method  have  remained,  more  or  less, 
unchanged  since  1988.  Therefore  we  selected  for  our  study  a  time  period  between  1988  and 
1993. 

During  the  period  selected  for  our  study  mining  works  were  carried  on  at  five  coal  levels 
simultaneously,  beneath  other  six  already  mined  out.  These  productive  levels  are: 

-  The  uppermost  productive  level  416  located  at  depths  of  700  to  740  meters.  Its 
thickness  ranges  from  4  to  5.5  m.  During  the  studied  period  mining  was  carried  on,  at  this  level, 
only  in  C  part  in  two  slabs; 

-  The  level  501  located  at  depth  range  from  770  to  810  m  with  thickness  ranging  from 
some  7  to  8.5m.  Coal  has  been  mined  there  in  three  slabs  in  NW,  C  and  NE  parts,  and  in  two 
slabs  in  SW  and  SE  parts,  where  the  coal  layer  was  thinner.  During  the  studied  period  mining 
was  carried  on  in  all  but  SE  parts.  Exploitation  of  level  501  was  the  most  endangered  by 
rockbursts.  Often  strong  tremors  caused  considerable  damages  to  mining  installations.  In  some 
cases  the  works  had  to  be  stopped  and  the  longwaUs  rebuild; 

-  The  level  504  located  from  2.5  to  15  meters  under  the  level  501.  Its  thickness  varies 
from  1.6  to  2.4  m.  The  level  in  its  productive  parts  was  nearly  mined  out  before  the  period 
selected  for  study.  Only  six  longwalls  were  working  in  the  years  1988-1993,  three  in  NW  part 
and  three  in  C  part; 

-  The  level  507.  It  appears  as  a  separate  coal  layer,  with  thickness  some  1.5  meters,  only 
in  NW  and  SW  parts  of  the  mine,  at  depth  of  some  745  m.  In  other  parts  it  joins  the  next  level 
510.  During  the  studied  period  NW  part  of  this  level  was  being  worked; 
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-  The  last  level  510  is  present  in  all  parts  of  the  mine  at  depth  range  from  800  to  850  m. 
The  coal  layer  is  thick,  5.3  to  6  m,  and  the  coal  has  been  mined  in  two  slabs.  In  the  years  1988- 
1993  the  works  were  done  in  NW,  SE  and  NE  parts  of  the  level. 

Up  to  50  meters  thick  layers  of  compact  rigid  sandstones  dominate  strata  separating  and 
overlaying  the  productive  layers.  The  other  strata  are  formed  of  shales  of  various  kind  with 
incidental  intrusions  of  coal. 

Maps  of  all  slabs  of  all  productive  levels  are  given  in  Appendices  1-9.  Dark  areas  denote 
parts  worked  out  before  1988.  Geometry  of  slopes  active  during  the  period  selected  for  study, 
marked  by  thick  black  lines,  is  combined  with  longwaU  positions  in  selected  days. 
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3.  Seismological  database 


A  local  mine  area  wide  seismic  network  was  installed  in  1975.  Since  1982  the  mine  has 
been  equipped  with  a  12  channel  distal  recording  system  called  PCM  G3,  produced  in  Warsaw, 
Poland.  The  sensors  (only  vertical  component)  have  been  distributed  all  over  the  mine  area  at 
various  depths  ranging  from  3(X)  to  850m  below  the  surface.  Horizontal  distances  between  the 
sensors  do  not  exceed  3.5  km.  An  overall  dynamics  of  the  recording  system  is  some  70  dB 
while  the  frequency  range  of  flat  response  is  form  0. 1  to  27  Hz. 

The  network  records  events  of  energy  from  lO^J.  The  recorded  signals  are  immediately 
processed  to  evaluate  their  parameters  and  then  stored  in  the  mine  catalog.  Every  event  in  the 
catalog  is  described  by  a  time  of  its  occurrence,  coordinates  of  its  epicenter,  its  seismic  soiffce 
energy,  comments  concerning  the  stope  and  the  level  where  the  event  occurred  and  additional 
comments  on  mining  works  to  that  could  have  an  influence  on  event  creation  (productive 
blasting,  destressing  blasting,  water  injection  etc.).  Various  test  have  pointed  out  that  the  average 
accuracy  of  epicenter  location  is  some  50  m  while  the  precision  of  energy  determination  is 
roughly  the  half  of  order.  Due  to  a  flat  distribution  of  sensors,  the  source  depth  cannot  be 
estimated  rehably.  The  only  information  about  vertical  positions  of  sources  comes  from  the 
comments  linking  events  to  i^ecific  coal  levels.  In  this  connection  we  consider  only  vertical 
division  of  the  data  into  different  productive  levels. 

For  some  of  the  weakest  events,  usually  of  energy  below  lO^J,  even  the  epicenter 
coordinates  are  not  determined.  We  use  these  data  only  in  descriptive  statistics  analyses.  For 
more  sophisticated  studies  we  need  values  of  all  parameters  of  events. 

14297  events  were  recorded  during  the  period  selected  for  the  study.  Out  of  them  1524 
occurred  at  the  level  416,  7365  at  the  level  501,  642  at  the  level  504,  1345  at  the  level  1345, 
2980  at  the  level  510,  and  for  441  the  source  location  was  not  determined.  The  distribution  of 
events  among  the  mine  parts  is  given  in  Figure  3.1.  Zero  level  denotes  events  for  which 
information  about  coal  level  was  not  provided.  SW  part  of  the  level  501  was  the  most 
seismically  active.  Also  E  part  at  the  level  510,  C  part  at  the  level  501,  C  part  at  the  level  416 
and  NW  part  at  the  level  507  experienced  considerable  activity. 


DISTRIBUTION  OF  SEISMIC  ACTIVITY  IN  WUJEK  MINE 
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As  mentioned  previously  the  catalog  comprised  events  of  energy  from  lO^J.  The 
strongest  recorded  event  had  energy  of  5.4*10*J  (Richter’s  magnitude  some  3.6).  It  occurred  at 
the  coal  level  501,  by  Arkona  fault,  well  apart  form  mining  works  of  longwall  I.  The  cumulative 
frequency-  energy  distribution  of  all  events  is  presented  in  Figure  3.2.  The  distribution  is  more 
or  less  linear  from  10  J.  Hence,  as  far  as  energy  is  concerned,  we  can  assume  that  the  database 
is  complete  starting  from  this  energy  value.  No  significant  differences  in  the  distributions 
constructed  for  activity  in  each  year  can  be  seen.  The  fact  suggests  a  time  uniformity  of  energy 
evaluation. 

The  spatial  distributions  of  seismicity  associated  with  the  level  416  during  the  studied 
period  together,  with  contours  of  active  excavations,  are  presented  in  Figures  3.3,  3.4,  3.5. 
Respective  figures  show  positions  of  sources  of  all  recorded  events,  of  events  that  energy  was 
greater  than  or  equal  to  10“*!  and  of  events  that  energy  was  greater  than  or  equal  to  lO^J.  Most  of 
events  occurred  close  to  active  stopes.  There  are,  however,  two  opposite  cases.  The  cluster 
located  in  SW  part  could  be  eventually  a  far  effect  of  severe  activity  at  the  501  level  .  .Similar 
phenomenon  at  the  same  place  is  observed  in  the  seismicity  map  of  the  levels  504  (Fig.  3.9)  and 
507  (Fig.  3.12).  Since  in  this  area  both  levels  (416  and  501)  are  distanced  vertically  of  some  60 
meters  the  activity  at  the  level  416  could  be  an  example  of  a  chaotic  dynamic  system  behavior  of 
rockmass  in  the  subcritical  state.  Such  hypotheses  explaining  strange  appearance  of  both  natural 
as  well  as  induced  seismicity  can  be  found  in  a  number  of  references  (  e.g.,  Grasso  1993, 
Morrison  1993).  The  other  cluster  located  in  NW  part  cannot  be  linked  neither  to  active 
excavations  nor  to  other  seismically  active  areas  at  any  productive  level. 

Generally  the  activity  of  level  416  during  the  studied  period  was  not  very  hazardous. 
There  were  only  three  events  of  energy  above  lO'*!  order. 

The  maps  of  seismic  sources  recorded  at  the  level  501  are  shown  in  Figures  3.6,  3.7,  3.8 
respectively.  The  most  active  areas,  with  many  strong  tremors,  were  connected  with  mining  in 
third  slab  of  SW  part,  adjoining  Arkona  fault,  and  with  mining  in  C  part.  The  cluster  whose 
center  is  at  some  (18600,  -8600)  cannot  be  linked  to  any  mining  or  seismic  activity  at  any  level. 

Figures  3.9,  3.10,  3.11  show  maps  of  seismic  sources  recorded  at  the  level  504.  Apart 
from  clusters  linked  to  excavations  in  C  and  NW  parts  one  may  see  a  group  of  seismic  events  in 
SW  part  of  possible  similar  origin  to  those  at  the  level  416  in  this  part.  This  seismic  activity  is. 
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Figure  3.4  Spatial  distribution  of  seismic  events  of  energy  >  10'*  J  associated  with  the  coal 
level  416 
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Figure  3.5  Spatial  distribution  of  seismic  events  of  energy  >  10*  J  associated  with  the  coal 
level  416 
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Figure  3.7  Spatial  distribution  of  seismic  events  of  energy  >  10‘‘  J  associated  with  the  coal 
level  501 
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however,  much  weaker  than  that  at  the  level  416,  although  the  level  504  is  much  closer  to  the 
level  501  than  the  level  416.  Some  slight  activity  in  the  northern  part  of  C  area  also  cannot  be 
directly  explained  by  mining  works.  In  general  energy  of  events  was  not  very  high.  There  were 
five  events  of  energy  above  the  lO"*!  order,  all  close  to  excavations  in  C  part. 

The  maps  of  seismic  sources  recorded  at  the  level  507  are  shown  in  Figures  3.12,  3.13, 
3.14.  The  level  exists  only  in  the  NW  and  SW  parts  but,  to  facihtate  comparisons  with  other 
levels,  the  same  scale  of  graphs  has  been  maintained.  The  seismic  activity  concentrated  in  areas 
of  active  excavations.  Events  probably  induced  by  seismic  activity  at  the  level  501  are  seen  in 
SW  part.  It  shoiild  be  noted  that  the  level  507  is  distanced  form  the  level  501  of  nearly  60 
meters  in  this  part  of  the  mine.  Despite  considerable  in  number  the  overall  seismicity  associated 
with  the  level  507  was  not  very  dynamic.  Only  nine  events  of  energy  above  the  lO^J  order  were 
recorded. 

The  maps  of  spatial  distributions  of  seismic  events  at  the  last  level  510  are  shown  in 
Figures  3.15,  3.16,  3.17.  One  can  see  clusters  connected  with  excavations  in  NW,  SE  and  NE 
parts.  Many  events,  however,  were  located  beyond  the  active  stopes,  particularly  in  NE  part  of 
the  mine.  Some  part  of  them  could  be  correlated  with  mining  works  some  35  meters  above,  in 
the  second  slab  of  coal  level  501 .  Still  an  open  question  remain  why  this  area  was  more  active  at 
the  level  510  than  at  the  productive  level  501.  Also  clusters  with  centers  located  at  some  (18150, 
-9(K))  and  (17950,  -96(X))  cannot  be  directly  linked  to  any  excavations  or  any  other  seismically 
active  parts  of  the  mine.  The  level,  in  general,  was  the  next,  after  the  level  501,  most  active  in 
high  energy  range.  Many  events  of  energy  starting  from  lO^J  occurred  in  its  SE  part. 
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Figure  3.12  Spatial  distribution  of  all  seismic  events  associated  with  the  coal  level  507 
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Figure  3.17  Spatial  distribution  of  seismic  events  of  energy  ^  10^  J  associated  with  the  coal 
level  510 
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4.  Descriptive  statistics 


The  aim  of  the  descriptive  statistics  studies  is  the  initial  recognition  of  statistical 
properties  of  database.  We  intend  to  check  whether  the  mining  induced  seismicity  originated  in 
various  parts  of  the  mine  is  so  different  that  simple  statistical  tests  can  reveal  these  differences. 
Similar  procedures  ought  to  be  used,  in  our  opinion,  while  studying  every  large  catalog  of  data, 
to  divide  it  into  coarsely  homogenous  parts  before  more  developed  methods  are  applied. 
Otherwise  strong  differences  will  be  found  by  the  methods  designated  to  study  detailed  structure 
of  the  data  and  they  can  dominate  results  making  detailed  analyses  impossible.  Additionally  the 
descriptive  statistics  comparisons  of  catalog  parts  firom  different  time  periods  give  chance  to  test 
indirectly  completeness  and  an  instrumental  time  homogeneity  of  the  catalog.  If  we  observed 
systematic  time-ordered  differences  in  distributions  of  a  parameter  (e.g.  consecutively  increasing 
mean  value)  we  could  suspect  diat  the  data  acquisition  changed  over  time. 

Since  major  features  of  spatial  distributions  of  sources  are  clearly  visible  from  the  maps 
(Figs.  3-17)  here  we  concentrated  ourselves  onto  time  of  evait  occurrence  and  event  energy. 
The  latter  parameter  is  strongly  non-linear.  TechnicaUy  it  is  much  easier  to  study  its  logarithm. 
Conclusions  drawn  for  the  logarithm  of  energy  obviously  are  also  valid  for  energy. 

As  fer  as  time  of  evait  occurrence  is  concerned  we  studied  the  interoccurrence  time.  For 
the  pure  Poissonian  process  of  event  generation  the  interoccurrence  time  has  exponential 
distribution.  It  is,  however,  known  that  inducing  processes  in  mines  are  not  purely  Poissonian 
(Lasocki  1992)  and  hence  empirical  distributions  of  this  parameter  does  not  match  the  theoretical 
exponential  distribution.  The  feature,  that  does  not  agree  with  main  event  series  behavior  in 
earthquake  seismology  (Gardner  and  Knopoff  1974),  comes  from  controlling  influence  of  time- 
varying  mining  works  over  mining  induced  processes.  The  non-Poissonian  structure  of  the 
induced  seismicity  was  ascertained  in  full  seismicity  series.  When  the  series  were  reduced  to 
only  strong  events  (E>10*J)  it  became  Poissonian  (Lasocki  1992). 

Our  descriptive  statistics  studies  mainly  concern  testing  differences  in  central  tendency 
and  empirical  distributions  of  parameters.  Since  distributions  of  both  the  interoccurrence  time  as 
well  as  the  logarithm  of  energy  are  unknown  but  surely  are  not  normal  we  basically  use 
nonparametric  testing  procedures.  These  are: 
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For  the  central  tendency: 

-  The  Mann-Withney  U  test.  The  test  evaluates  differences  in  medians  between  two 
samples.  The  test  is  computed  based  on  rank  sums  of  the  data.  It  is  regarded  as  the  most 
powerful  nonparametric  alternative  to  the  t-test  for  testing  differences  between  means  of 
normally  distributed  samples.  The  test  was  even  found  to  be  more  powerful  than  the  t-test  in 
some  cases; 

The  KruskaU-Wallis  ANOVA  by  Ranks  and  the  median  tests.  The  procedures  are 
nonparametric  alternatives  to  <Hie-way  analysis  of  variance  and  test  for  significant  differences  in 
medians.  Having  more  than  two  samples  to  compare  the  use  of  these  tests  may  be  preferred  over 
the  univariate  approach; 

For  the  distribution; 

-  The  Kolmogorov-Smimov  two  sample  test  that  assesses  the  hypothesis  that  two  samples 
were  drawn  from  different  populations.  Unlike  the  tests  for  central  t^dency  it  is  also  sensitive  to 
differences  in  the  general  shapes  of  the  distributions  in  the  two  samples; 

-  The  Wald-Wolfowitz  runs  test.  The  test  assesses  the  hypothesis  that  two  samples  were 
drawn  from  two  populations  that  differ  with  respect  to  the  central  tendency  and  to  the  general 
shape  of  the  distribution. 

The  detailed  description  of  these  tests  can  be  found  in  a  number  of  books  on  mathematical 
statistics  (e.g.,  Hollander  and  Wolfe  1973,  Box  et  all  1987). 

Having  in  mind  remarkable  robustness  to  deviations  from  normality  of  the  standard 
analysis  of  variance  sometimes,  to  cross-vahdate  our  conclusions  implied  by  nonpiarametric 
procedures,  we  use  the  one-way  analysis  of  variance  to  study  differences  in  means.  We  maVp 
use  of  the  LSD  (Least  Significant  Differences)  test. 

All  mentioned  procedures  provide  the  probability  of  error  when  accepting  the  hypothesis 
about  the  existence  of  the  difference.  The  probability  is  called  the  significance  level.  One  usually 
accepts  the  hypothesis  when  the  significance  level  is  less  than  0.05.  However  when  the 
probability  is  greater  than  this  value  one  cannot  conclude  that  the  samples  were  drawn  Ifom  the 
same  population  (nor  that  they  have  the  same  central  tendencies).  If  the  probability  of  error 
when  accepting  the  hypothesis  about  the  existence  of  the  difference  is  say  0.1  the  probability  of 
error  when  accepting  the  hypothesis  about  equality  is  0.9. 
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Studies  of  the  logarithm  of  seismic  energy  of  events 

Differences  in  location  parameters  (mean,  median)  of  distributions  of  logarithm  of 
seismic  energy  between  samples  containing  events  recorded  at  different  coal  levels  in  particular 
years  of  the  studied  period  were  studied  by  means  of  the  Kruskal-Wallis  test  and  the  LSD  test. 
The  Kruskal-Wallis  test  proves  significance  of  the  differences  between  parts  of  the  seismicity 
series  due  to  the  considered  years,  for  all  seismicity  recorded  in  the  whole  mine,  as  well  as  for 
seismicity  connected  with  the  specific  coal  levels.  The  differences  are  veiy  highly  significant. 
The  level  of  significance  for  all  studied  groups  is  less  than  10"*. 

The  detailed  insight  into  these  differences  and  the  differences  in  the  general  shape  of  the 
logE  distributions  is  assessed  with  the  use  of  the  LSD  test  for  means  and  the  Kolmogorov- 
Smimov  two  sample  test,  req)ectively.  The  logE  distribution  is  not  normal  but  violation  of  the 
normality  assumption  has  no  significance  for  the  LSD  test  results  when  samples  are  large.  Table 
4.1  presents  results  of  the  two  tests.  The  results  of  the  Kolmogorov-Smimov  testing  procedure 
has  been  fully  confirmed  by  the  Wald-WoUbwitz  test. 

Similarities  (non-significant  differences)  in  both  mean  value  and  shape  are  randomly 
distributed  among  years  of  the  studied  period,  however  the  slight  tendency  to  decrease  logE 
mean  value  with  time  can  be  picked.  In  some  cases  the  samples  that  have  similar  mean  have  also 
similar  shape  of  the  energy  distribution,  in  other  they  do  not. 

The  differences  in  both  the  shape  and  the  location  of  energy  distributions  between 
seismic  series  recorded  in  various  years  prove  a  strong  link  between  seismic  energy  generation 
and  mining  works.  The  works  changed  in  time  both  their  location  and  their  intensity.  Some 
stopes  were  closed,  in  some  other  parts  the  exploitation  was  commenced.  Irregular  variations  of 
mining  activity  caused  irregularities  in  seismicity  generation.  Valuable  information  about  a 
nature  of  mining  works  influence  onto  seismicity  generation  can  be  extracted,  in  our  opinion, 
after  confronting  the  obtained  results  with  detailed  data  about  mining  works  at  particular  levels. 
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Table  4.1.  Significance  of  differences  between  means  and  distributions  of  logE  for  samples  from  different  years 


X  LSD  test 


Year 

Coal  level 

Sample  size 

Mean  value 

Significance  level 

1988 

1989 

1990 

1991 

1988 

all  levels 

4038 

3.69 

X 

0.13 

0.03 

0.98 

416 

1269 

3.57 

X 

0.003 

0.54 

- 

501 

1670 

3.78 

X 

0.0001 

<10-® 

0.(K)07 

504 

239 

3.55 

X 

0.001 

0.004 

0.11 

507 

525 

3,71 

X 

2x10'^ 

0.07 

0.0002 

510 

194 

3.74 

X 

0.99 

0.005 

0.04 

1989 

all  levels 

3219 

3.70 

>0.050 

X 

0.0009 

0.23 

416 

219 

3.45 

<0.001 

X 

0.08 

- 

501 

2024 

3.71 

<0.001 

X 

0.018 

1.0 

504 

98 

3.34 

<0.050 

X 

10“® 

0.0006 

507 

438 

3.85 

<0.001 

X 

<10-® 

<10^ 

510 

326 

3.74 

>0.050 

X 

0.001 

0.02 

1990 

all  levels 

2052 

3.65 

<0.001 

<0.001 

X 

0.08 

416 

34 

3.63 

>0.050 

<0.050 

X 

- 

501 

1380 

3.66 

<0.001 

<0.001 

X 

0.03 

504 

73 

3.75 

<0.001 

<0.001 

X 

0.61 

507 

224 

3.63 

<0.050 

<0.001 

X 

0.022 

510 

299 

3.60 

<0.001 

<0.010 

X 

0.58 

1991 

all  levels 

1680 

3.69 

<0.001 

<0.001 

<0.001 

X 

416 

0 

0 

- 

- 

- 

X 

501 

1297 

3.71 

<0.001 

<0.001 

<0.001 

X 

504 

38 

3.69 

<0.001 

<0.001 

>0.050 

X 

507 

86 

3.48 

<0.001 

<0.001 

<0.050 

X 

510 

204 

3.63 

<0.001 

<0.010 

>0.050 

X 

1992 

all  levels 

1371 

3.60 

<0.001 

<0.001 

<0.001 

<0.001 

416 

0 

0 

- 

- 

- 

- 

501 

453 

3.55 

<0.001 

<0.001 

<0.010 

<0.001 

504 

101 

3.58 

>0.050 

<0.010 

<0.050 

<0.050 

507 

39 

3.49 

>0.050 

<0.001 

>0.050 

>0.050 

510 

703 

3.62 

<0.010 

<0.001 

<0.010 

>0.050 

1993 

aU  levels 

1940 

3.58 

<0.001 

<0,001 

<0.001 

<0.001 

416 

0 

0 

- 

- 

- 

- 

501 

527 

3.64 

<0.001 

<0.001 

<0.001 

<0.001 

504 

92 

3.54 

<0.010 

>0.050 

<0.050 

>0.050 

507 

26 

3.49 

>0.050 

<0.010 

>0.050 

>0.050 

510 

1252 

3.55 

<0.001 

<0.001 

<0.001 

<0.050 

Kolmogorov-Smimov  test 

Qigl] 
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The  general  trend  to  lower  mean  energy  in  time  could  be  due  to  a  general  decrease  of 
number  of  strong  events  in  time  during  the  studied  period.  However,  it  is  also  possible,  that  the 
decreasing  trend  could  be  resulted  by  an  enhancement  of  recording  facilities  in  the  mine.  In  this 
case,  increasing  in  time  number  of  weak  events  recorded  could  reduce  the  mean  value.  If  the 
latter  was  true  then  the  instrumental  time-homogeneity  of  the  database  would  not  be  preserved. 
The  problem  can  be  resolved  through  interoccurrence  time  studies.  If  the  latter  explanation  is 
valid  then  a  tendency  to  decrease  interoccurrence  time  (increase  the  mean  event  rate)  in 
consecutive  years  will  be  noted. 

In  order  to  study  the  same  as  previously  properties  of  logE  distributions  in  the  series 
comprised  (mly  strong  events  we  removed  from  the  samples  all  events  weaker  than  3.2x10"^  J 
(half  of  order  10^  and  repeated  the  analysis  of  differences  between  parameters  of  location.  The 
sample  reduction  decreases  considerably  sample  sizes.  In  this  case,  due  to  not  normal 
distribution  of  logE,  only  nonparametric  procedures  can  be  used.  The  a^tlied  Krudcall-Wallis 
and  the  median  tests  reject  hypothesis  about  significant  differences  between  the  reduced  samples 
associated  with  particular  years  of  the  studied  period.  The  result  pointes  out  that  the  strong  event 
series  are  less  variable  in  time,  concerning  energy  of  events,  than  the  full  series.  Since  we 
suspect  that  the  variation  of  energy  distribution  in  time  is  due  to  the  variation  of  mining  works 
we  can  conclude  that  generation  of  seismicity  in  higher  energy  range  is  less  dependent  upon  the 
mining  activity.  This  gives  a  strong  suggestion  that  a  control!  of  mining  works  over  at  least  a 
considerable  part  of  strong  event  generation  is  weaker.  Thus  the  other  than  directly  controlled  by 
mining  process  takes  place.  As  the  things  are  an  important  conclusion  may  be  put  forward: 
There  are  probably  more  than  one  inducing  processes  and  our  induced  seismicity  series 
comprises  mixture  of  their  outcomes. 

Besides  time  divided  series  studies  also  series  due  to  particular  coal  levels  and  due  to 
praticular  areas  of  the  mine  are  analysed.  Because  of  limited  sample  sizes  we  use  here  tiie  Mann- 
Withney  U  test  to  compare  medians.  The  results  of  testing  differences  between  samples 
associated  with  the  coal  level  in  the  respective  mine  areas  are  presented  in  Table  4.2,  whereas 
those  between  samples  associated  with  the  mine  area  at  the  respective  coal  levels  are  shown  m 
Table  4.3.  The  results  for  full  samples  are  above  the  diagcmals  of  the  tables  while  the  values 
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assemblyed  below  the  diagonals  concern  samples  with  events  of  energy  equal  to  or  greater  than 
lO'^J.  Relations  only  between  samples  of  satisfactoiy  size  are  presented. 

As  can  be  seen  in  the  Table  4.2  the  locations  (medians)  of  logE  distributions  do  not 
differ  significantly  between  seismic  series  at  the  levels  416  and  504  and  between  the  series  at  the 
levels  501  and  510.  The  first  pair  of  levels  comprises  thin  almost  exploited  coal  layers. 
Seismicity  associated  with  them  neither  had  big  rate  nor  was  very  dynamic.  Coal  layers  of  the 
second  pair  are  thick  and  mining  there  was  accompanied  by  considerable  seismic  activity  in  wide 
energy  range. 

The  strong  event  series  studies  form  another  pairs  with  non-diffeent  logE  locations. 
These  are:  (416,  507) ,  (501,  504)  and  (504,  510).  The  links  caimot  be  explained  on  the  basis  of 
layer  characteristics.  The  strong  event  generation  is  probably  less  dependent  upon  coal  layers 
features  but  more  controlled  by  local  stress  field  due  to  structures  located  around  the  stopes  and 
in,  above  and  below  the  coal  seam. 

In  general  differences  between  distribution  locations  of  the  considered  mine  areas  are 
stronger  in  full  samples  than  in  reduced  data.  The  full  seismic  series  from  NW  part,  at  the  level 
507  differs  firom  all  other  full  series  from  diis  part,  at  all  other  levels.  Significant  differences  are 
also  found  between  the  full  series  from  the  levels  504  and  510  in  NW  part.  Only  two  levels:  504 
and  510  in  C  part  are  not  different  from  the  point  of  view  of  full  series  logE  distributions. 

The  differences  decreases  after  the  energy  threshold  is  applied.  Neither  full  nor  reduced 
series  recorde  in  SW,  NE  and  SE  parts  at  different  coal  levels  differ  significantly. 
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Table  4.2.  Significance  of  differences  between  medians  of  logE  for  samples  from  different  coal  levels 


X  FULL  SAMPLE  SIGNIFICANCE  LEVEL 


Coal  level 

Region 

Mean  in 
full  sample 

Mean  in 
reduced  sample 

416 

501 

504 

507 

510 

416 

whole  level 

3.56 

4.25 

X 

<10-* 

0.42 

<10*^ 

<10-^ 

SW 

3.58 

4.24 

X 

0.14 

- 

0.07 

0.14 

C 

3.58 

4.25 

X 

<10-* 

0.04 

- 

0.046 

501 

whole  level 

3.70 

4.42 

<10^ 

X 

<10*^ 

7x10'^ 

0.37 

NW 

3.62 

4.25 

- 

X 

0.78 

0.0002 

0.22 

SW 

3.77 

4.43 

0.05 

X 

- 

0.35 

0.89 

c 

3.81 

4.41 

3x10^ 

X 

<10-® 

- 

8x10^ 

NE 

3.78 

4.25 

- 

X 

- 

- 

0.5 

SE 

3.77 

- 

- 

X 

0.005 

- 

0.84 

504 

whole  level 

3.55 

4.39 

0.0008 

0.83 

X 

o 

V 

5x10’^ 

NW 

3.59 

4.25 

- 

0.60 

X 

<10-® 

0.03 

C 

3.66 

4.41 

0.0002 

0.52 

X 

0.96 

SE 

3.48 

- 

- 

- 

X 

- 

0.0003 

507 

whole  level 

3.72 

4.25 

0.78 

o 

V 

0.0004 

x 

<10*^ 

NW 

3.77 

4.25 

- 

0.36 

0.97 

X 

0.0001 

SW 

3.78 

4.31 

0.17 

0.64 

- 

X 

0.6 

510 

whole  level 

3.61 

4.34 

0.01 

7x10'^ 

0.06 

0.005 

X 

NW 

3.69 

4.33 

- 

0.98 

0.77 

0.26 

X 

C 

3.64 

4.30 

0.35 

0.12 

- 

X 

NE 

3.74 

4.23 

0.82 

_ 

X 

REDUCED  SAMPLE  SIGNIFICANCE  LEVEL 
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Full  seismic  series  from  various  parts  of  the  mine  (Table  4.3)  turn  out  to  be  similar  (not 
different)  in  a  difficult  to  explain  way.  More  obvious  similarities,  that  can  be  correlated  with 
geographical  distribution  of  the  studied  areas  are  found  for  reduced  samples.  These  similar  pairs, 
from  the  point  of  view  of  seismic  energy  generation,  are:  (SW,  C) ,  (NW,  NE)  and  (SE,  Q.  It 
is  then  likely  that  the  full  series  generation,  in  which  weak  events  dominate,  is  more  level 
controlled  while  the  stong  event  series  generation  is  more  region  controlled.  At  particular  levels 
the  same  blurring  of  differences  between  the  seismic  series  from  different  parts  of  the  mine  after 
applying  energy  threshold  can  be  observed  in  Table  4.3  as  it  was  noted  for  the  series  from 
different  coal  levels  at  particular  regions  in  table  4.2. 


Table  4.3.  Significaiice  differences  between  medians  td'lo^  for  samples  from  different  mine  parts  (regions) 


X  FULL  SAMPLE  SIGNIFICANCE  LEVEL 


Region 

Coal  level 

Mean  in 
liill  sample 

Mean  in 
reduced  sample 

NW 

SW 

c 

NE 

SE 

NW 

all  levels 

3,55 

4.26 

X 

0.56 

<10-^ 

0,043 

0.18 

501 

3.62 

4.25 

X 

0.011 

0.002 

0.001 

0.11 

504 

3.59 

4.25 

X 

- 

0.38 

- 

0.02 

507 

3.77 

4.25 

X 

0.91 

- 

- 

- 

510 

3.69 

4.33 

X 

0.65 

0.27 

0,004 

0.12 

SW 

all  levels 

3.77 

4.43 

A 

o 

X 

<10-® 

0.011 

0.27 

416 

3.68 

4,24 

- 

X 

0.014 

- 

- 

501 

3.77 

4,43 

0.38 

X 

0.0007 

0,13 

0.89 

507 

3.78 

4.31 

0.09 

X 

- 

- 

- 

510 

3.74 

- 

- 

X 

0.39 

0.57 

0,99 

C 

all  levels 

3.71 

4.38 

<10-^ 

0.14 

X 

<10"^ 

0.01 

416 

3.58 

4,25 

0.60 

X 

- 

- 

501 

3.81 

4.41 

0.35 

0.85 

X 

0.82 

0.56 

504 

3.66 

4,41 

0.14 

- 

X 

- 

0.015 

510 

3.65 

4.30 

0.94 

- 

X 

0.0002 

0.008 

NE 

all  levels 

3.75 

4.23 

0.75 

5x10'^ 

7x10'^ 

X 

0.006 

501 

3.79 

4.25 

0.49 

0.08 

0.06 

X 

0.41 

510 

3.74 

4.25 

0.44 

- 

0.55 

X 

0.57 

SE 

all  levels 

3.73 

4.37 

0.003 

0.03 

0.28 

0.02 

X 

510 

3.73 

4.37 

0.40 

0.57 

0.072 

X 

REDUCED  SAMPLE  SIGNIFICANCE  LEVEL 
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Interoccurrence  time  studies 


Results  of  the  Mann-Whitney  U  test  of  differences  between  interoccurrence  time  medians 
of  events  recorded  in  the  whole  mine  area  in  particular  years  of  the  studied  period  are  presented 
in  Table  4.4.  In  all  but  one  case  differences  are  very  highly  significant.  The  medians  show  clear 
tendency  of  incrasing  mean  interoccurrence  time  during  the  period  from  1988  to  1992.  The 
central  tendencies  of  the  interoccurence  times  distributions  represent  the  reciprocal  of  central 
tendencies  of  the  event  rates.  Thus  we  reject  the  possibility,  discussed  in  the  previous  subsection, 
that  the  seismicity  time  variation  could  be  due  to  raising  recording  ability  of  the  seismic  network 
in  the  studied  period.  If  this  was  true  we  would  observe  a  reversed  teidency  of  the 
interocurrence  time  median  caused  by  a  general  increase  of  number  of  recorded  events.  We  can 
then  attribute  both  the  logE  distributions  and  the  interoccurrence  times  distributions  tendencies  to 
the  decrease  of  overal  activity  during  the  studied  period.. 


Table  4.4.  Results  of  median  test  for  interoccurrence  time  of  events  recorded  in  the  Wujek  mine 


Year  Sample  size  Median  Signiflcance  level 


1989 

1990 

1991 

1992 

1993 

1988 

4038 

77 

<10-® 

<10-® 

<10-® 

<10'® 

<10® 

1989 

3219 

97 

X 

<10-® 

<10-® 

<10‘® 

<10® 

1990 

2052 

134 

X 

<10-® 

<10-® 

0.46 

1991 

1680 

176 

X 

0.0001 

<10® 

1992 

1371 

208 

X 

<10® 

1993 

1940 

125 

X 

Results  of  the  U  test  of  differences  between  interoccurrence  time  medians  for  events 
recorded  at  the  particular  mining  levels,  in  particular  years  of  the  studied  period  are  shown  in 
Tables  4.5, 4.6, 4.7, 4.8,  4.9,  respectively. 


Table  4.5.  Results  of  median  test  for  interoccurrence  time  of  events  recorded  in  the  coal  level  416 


Year 

Sample  size 

Median 

Significance  level 

1989 

1990 

1988 

1272 

145 

<10-^ 

<10-® 

1989 

219 

969 

X 

0.0001 

1990 

34 

3025 

X 

Table 

4.6.  Results  of  median  test  for  interoccurrence  time  of  events  recorded  in  the  coal  level  501 

Year 

Sample  size 

Median 

Significance  level 

1989 

1990 

1991 

1992 

1993 

1988 

1673 

1% 

<10-® 

0.012 

3x10-® 

<10-® 

<10^ 

1989 

2027 

154 

X 

<10^ 

<10-® 

<10-® 

<10^ 

1990 

1384 

215 

X 

0.04 

<10-® 

<10-® 

1991 

1299 

250 

X 

0.0001 

<10-® 

1992 

454 

762 

X 

6x10® 

1993 

529 

528 

X 

Table  4.7.  Results  of  median  test  for  interoccurrence  time  of  events  recorded  in  the  coal  level  504 

Year  Sample  size 

Median 

Significance  level 

1989 

1990 

1991 

1992 

1993 

1988 

239 

378 

o 

V 

10-® 

<10^ 

<10-® 

<10-® 

1989 

98 

1792 

X 

0.44 

0.36 

0.68 

0.56 

1990 

73 

1393 

X 

0.17 

0.57 

0.17 

1991 

38 

2306 

X 

0.33 

0.82 

1992 

101 

1990 

X 

0.32 

1993 

93 

2151 

X 

Table  4.8.  Results  of  median  test  for  interoccurrence  time  of  events  recorded  in  the  coal  level  507 


Year 

Sample  size 

Median 

Significance  level 

1^ 

1990 

1991 

1992 

1993 

1988 

528 

538 

0.40 

0.12 

<10^ 

<10-® 

<10‘® 

1989 

440 

608 

X 

0.34 

<10-* 

A 

O 

O 

V 

1990 

227 

661 

X 

<10-* 

.3x10-^ 

<10"® 

1991 

86 

2253 

X 

0.34 

0.09 

1992 

39 

4338 

X 

0.68 

1993 

26 

4772 

X 

Table  4.9.  Results  of  median  test  for  interoccurrence  time  of  events  recorded  in  the  coal  level  510 


Y'ear  Sample  size  Median  Signiflcance  level 


1989 

1990 

1991 

1992 

1993 

1988 

194 

886 

0.04 

0.09 

0.14 

o 

V 

o 

V 

1989 

326 

799 

X 

0.72 

8x10'* 

<10-* 

<10-* 

1990 

299 

860 

X 

0.0004 

<10-® 

<10‘® 

1991 

204 

1150 

X 

<10^ 

<10-^ 

1992 

703 

392 

X 

O 

V 

1993 

1254 

166 

X 

Significant  or  higMy  significant  differences  are  found  between  medians  of  interoccurence 
times  of  seismic  series  recorded  at  the  levels  416  and  501.  The  seismic  activity  of  the  level  501 
varied  irregularly. 

No  significant  differences  in  the  medians  occurres  between  the  series  from  the  level  504 
except  the  median  for  the  series  from  1988  that  was  considerable  lower  than  those  for  the  series 
from  other  years. 

Two  homogeneous  groups  of  series  concerning  locations  of  interoccorrence  times 
distributions  are  found  for  the  level  507.  One  comprises  the  seismic  series  recorded  in  the  years 
1988-1990  while  the  other  the  series  recorded  in  the  years  1991-1993  when  the  seismic  activity 
was  very  low. 

The  central  location  of  interoccurrence  time  distribution  of  events  associated  with  the 
level  510  tumes  out  to  be  not  different  in  the  years  1988-1990.  Exploitation  of  the  level  510  in 
slab  1  in  NE  part,  that  started  by  the  end  of  1992,  caused  considerable  increase  of  the  seismic 
activity  in  the  level. 


Conclusions 


Although  the  descriptive  statistics  studies  are  a  complementaiy  tool  in  our  work  they 
justify  the  following  suggestions  and  conclusions: 

1 .  Seismicity  is  a  complex,  both  time  and  place  dependent  phenomenon. 

2.  The  event  rate  and  the  energy  distribution  of  seismicity  is  time-varying  due  to  time  changes  of 
mining  works.  The  overal  process  tended  towards  lower  energy  and  lower  rate  during  the 
studied  period.  The  tests  do  not  evidence  instrumental  inhomegeneity  of  the  database. 

3.  The  distributions  of  energy  of  strong  events  do  not  vary  significantly  in  time.  This  proves  the 
existence  of  at  least  two  inducing  processes  in  the  mine.  The  one,  time-varying,  controlled  by 
mining  works,  generates,  in  general,  weaker  events.  The  other  one,  weakly  dependent  on  time, 
i.e.  mining  activity,  gives  rise  to  stronger  events. 

4.  Statistical  features  of  series  comprising  events  from  the  whole  energy  range  show  some 
correlation  with  prqrerties  of  coal  layers  they  are  associated  with.  Statistical  characteristics  of 
series  formed  by  only  strong  events  are  betta-  correlated  with  the  part  of  the  mine  they  came 
from  than  with  the  coal  level  they  were  associated  with.  The  process  of  generating  weak  events 
is  predominantly  dependent  upon  the  coal  layer  and  the  stope  characteristics.  The  process 
responsible  for  generating  stronger  events  is  more  dependent  upon  joint  influence  of  structures  in 
an  enviromnent  of  the  excavation  than  directly  dependent  upon  the  excavation  itself. 

5.  The  parameters  of  strong  event  series  from  various  mine  areas  are  less  different  tiian  the 
parameters  of  full  senes  from  these  areas.  The  generating  process  mainly  responsible  for  strong 
events  occurrences  is  less  place  dependent  compared  to  the  process  of  weak  event  generation. 
This  is  another  suggestiwi  that  the  process  of  strong  event  generation  is  controlled  by  some 
generlized  superposition  of  geologic  and  mining  structures  of  the  mine  area. 

6.  As  mentioned  the  descriptive  statistics  studies  enable  preliminaiy  data  recognition  thus  the 
conclusions  cannot  be  considered  as  fimte  and  fully  reliable.  In  particular,  we  are  not  allowed  to 
state  that  there  are  only  two  generating  processes.  Howevo",  the  studies  evidence  the  feet  that  the 
induced  seismicity  is  a  multimodal  process  and  its  multimodality  is  reflected  in  the  distribution 
of  evait  energy. 


5.  Maximum  value  of  enra’gy  studies 


Studies  of  statistical  properties  of  the  largest  annual  magnitudes  to  assess  characteristics 
of  earthquake  series  have  been  j^esented  in  the  number  of  papers  (e.g.,  Nordquist  1945,  Epstein 
and  Lomnitz  1966,  Yegulalp  and  Kuo  1974,  Lomnitz  1974,  Kijko  and  Sellevoll  1981,  Campbell 
1982,  Kijko  1982,  Gan  and  Tung  1983).  This  approach  can  also  be  met  in  the  field  of  mining 
seismicity  (e.g.,  Kushnir  et  al.  1984,  Dessokey  1984,  McGair  1984,  Lasocki  1993,  Gibowicz 
and  Kijko  1994,  Kijko  and  Funk  1994).  The  popularity  of  the  extreme  value  methods  stems 
from  the  lowered  requirements  concerning  the  recording  system,  that  should  be  able  to  record 
the  strongest  events. 

Let  Xj  ,  i=l...,n  be  outcomes  of  the  certain  process,  obtained  during  the  certain  time 
penod,  say  ft,  t+ATJ.  Let  the  cumulative  distribution  function  of  the  random  variable  X  be  F(x). 
When  we  divide  the  time  period  into  m  equally  long  intervals  At  {AT  =mAt)  and  we  find 
yk=tnwc(Xj\  t+(k-l)At  <t+kAt),  k=l,..,m  we  diall  have  outcomes  of  another  process  of 
maximum  generatitMi  with  its  own  cumulative  distribution  fimction,  say  G(y).  Gumbel  (1962) 
has  proved  that  there  are  only  three  possible  forms  of  G(y)  when  n  tends  towards  infinity.  These 
asymptotic  extreme  value  distributions  are: 
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where:  0  >  0  and  k  >  0. 


The  exact  distribution  of  the  extreme  value  T  can  be  reached  only  when  the  underlying 
distribution  of  the  random  variable  X  is  known.  The  specific  solutions  (e.g.,  Lasocki  1993, 
Gibowicz  and  Kijko  1994,  Kijko  and  Funk  1994)  usually  assume  a  kind  of  stationarity  of  the 
event  series,  which  leads  to  the  Poisson  distribution  of  the  number  of  events  in  a  time  interval. 
Since  induced  seismicity  series,  from  such  long  time  periods  as  studied  in  this  work,  are  not 
Poissonian  (Lasocki  1992)  we  had  to  use  the  asymptotic  forms  of  the  distribution  of  maxima. 

The  asymptotic  distribution  of  extrema  Y  is  related  to  the  distribution  of  the  random 
variable  X  (Johnson  and  Kotz  1970).  The  variable  Thas  type  I  asymptotic  distributkai  if; 


lim  «|l  -  G(X,_„  - 1 + ]  )}  =  exp(-;;) 


n  ->  00 

where  X^ :  G(XJ=z, 
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the  type  n  distribution  : 
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X  — >  00 
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where  c>0  and  k>0. 


and  type  III  distribution  if: 

1-G(cx^  , 

l-G(x  +  c.)  (5-6) 

X  0 

The  condition  (5.4)  is  satisfied  for  unbounded  distributions  like  normal  and  exponential, 
the  condition  (5.5)  for  distributions  bounded  from  below  and  the  condition  (5.6)  for  those 
bounded  form  above.  The  energy  of  induced  events  has  either  the  unlimited  distribution  assigned 
by  the  Gutenberg-Richter  relation  or  a  limited  from  above  form  (Cosentino  et  al.  1977,  Lasocki 
1993,  1993a,  Kijko  and  Funk  1994).  Thus  the  asymptotic  distribution  of  logarithm  of  energy 
maxima  is  either  of  type  I  or  of  type  IIL  In  the  first  case  the  logarithmic  function  governs  the 
relation  between  the  logarithm  of  energy  maxima  and  the  logarithmized  empirical  cumulative 
distribution  function  of  logE  maxima.  In  the  second  case  this  relation  is  the  power  function 
asymptotically  approaching  a  limiting  value  of  logE.  Every  violation  of  the  expected  shapes  by 
the  empirical  c.d.f.  of  energy  maxima  can  be  interpreted  as  multimodality  of  the  underlying 
energy  distribution,  that  is  the  sample  contains  mixture  of  outcomes  of  different  generating 
processes.  First  attempts  to  use  the  extreme  value  mefliods  to  investigate  complexity  of  the 
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seismic  catalogs  from  mines  were  given  by  Kijko  et  al.  (1987).  The  other  applications  for  the 
regional  seismicity  of  the  Upper  Silesian  Coal  Basin,  Poland  were  presented  by  Idziak  et  al. 
(1991).  The  method  would  provide  distinct  results  when  the  seismicity  generation  processes  have 
different,  well  separated  energy  modes  and  the  event  rate  of  the  high  energy  mode  is  lower  than 
the  others.  The  main  drawback  of  the  method  is  the  use  of  the  cumulative  distribution  function, 
that  is  the  integral  frmction,  which  screens  local  effects  and  insufficiently  large  peaks  of  the 
jM-obability  density  fimction  are  not  detectable.  The  method  requires  large  size  samples  to  reach 
the  asymptotic  distribution  of  maxima. 

Identification  of  the  seismicity  generation  processes  by  the  extremum  value  statistics  studies 

Based  on  the  theoretical  development  given  above  a  number  of  logE  maximum  value 
studies  were  made.  Figure  5.1  shows  probability  distribution  of  the  2()-day  maximum  seismic 
energy  release  at  the  whole  mine  area,  in  the  studied  period.  The  size  of  our  database,  in 
average  130  events  per  20-day  interval,  allowed  to  compare  the  empirical  distributions  with  their 
asymptotic  forms.  The  shape  of  the  empirical  distribution  function  is  smooth.  One  can  clearly 
see,  however,  at  least  two  points  of  inflexion.  These  two  most  distinct  points  evidence  at  least 
dnee-mode  structure  of  the  seismic  energy  distribution  of  all  events  from  the  mine. 
Consequently  they  prove  existence  of  at  least  three  different  generating  processes,  each  of 
different  energy  mode.  The  low  energy  mode  seems  to  tend  asymptotically  to  the  limiting  value 
of  energy  at  some  3x10^  J.  The  medium  energy  mode  has  the  upper  limit  at  some  lO’  J.  The 
high  energy  mode  limit  exceeds  10*  J. 

The  extreme  value  studies  require  the  statistical  structure  of  the  sample  of  maxima. 
There  cannot  be  any  significant  systematic  time  variations  of  the  maximum  distribution.  The 
homogeneity  of  our  maximum  value  set  was  tested  by  the  Kruskal-Walhs  and  median  tests  of 
differences  in  locations  and  by  the  Kolmogorov-Smimov  test  of  differences  in  the  shapes  of 
distributions.  Both  tests  rejected  hypotheses  about  significant  differences  between  the  maximum 
energy  distributions  of  seismic  events  series  recorded  in  the  particular  years  of  the  studied 
period. 

As  was  proved  by  the  descriptive  statistics  studies,  the  seismicity  is  a  highly 
nonhomogeneous  structure.  The  distributions  of  event  parameters  are  different  in  different  parts 


of  the  mine.  The  distribution  of  energy  maxima  shown  in  Figure  5.1  represents  a  generalized 
form  obtained  while  summing  up  local  effects.  There  are  secondary  modes  in  the  distribution, 
which,  because  of  this  summing,  cannot  be  picked.  To  detail  our  studies  we  investigated  the 
seismic  series  associated  with  the  particular  coal  levels  and  the  particular  mine  parts.  In  all 
further  discussed  cases  we  used  the  20-day  intervals  to  evaluate  the  logE  maxima.  The  number 
of  events  per  interval  was  set  to  be  greater  than  30,  hence  the  distributions  could  be  regarded  as 
approaching  their  asymptotic  forms.*  The  logE  maxima  from  the  intervals  for  which  this 
condition  was  not  satisfied,  were  not  used  to  build  the  empirical  distributions. 

The  probability  distributions  of  maximum  seismic  energy  release  at  the  particular  levels, 
in  the  studied  period,  are  shown  in  Figures  5.2  -  5.6,  respectively.  The  following  modes  of  the 
distributions  were  identified: 

-  level  416:  I  -  En^5xl0'*  J  ,  II  -  J  and  possibly  the  third  one  with  only  one 

observation; 

-level  501:  possibly  I  -  E^5xl0^  J  ,  n-E^^lO**!  ,  m-E^5xlO‘*J  , 
IV-E^>10*J; 

-  level  504:  I  -  ETOax*’5xl0'*  J  ,  II  -  E„,„,^10^  J.  Limited  number  of  observations 

makes  this  identification  uncertain; 

-level  507:  possibly  I  -  £^5x10^*  J  ,  n-E^«2xlO*J; 

-level  510:  possibly  I  -  E^W^SxlO'*  J  ,  n-E^«2xlO^J,  HI  -  E„,3„«5xl0®  J  and 
possibly  IV  -  Ejnax>  10*  J  with  only  one  observation. 

The  mode  with  limiting  value  in  the  first  half  of  order  10^  J  is  present  in  distributions  of  logE 
maxima  for  seismic  series  from  every  level,  although  for  the  distribution  associated  with  the 
level  501  it  is  obscured  by  the  higher  energy  modes.  Some  traces  of  the  possible  mode  of  logE 
maximum  distribution,  with  limiting  value  «5xl0‘*  J  is  also  noticed  for  the  series  from  every 
level  except  for  that  from  the  level  501.  The  levels  501  and  510  associated  distributions  display 


'  When  the  process  of  event  occurrence  is  Poissonian,  the  exact  distribution  of  maxima  can  be  rqrlaced  by  its 
asymptotic  form  if  exp(-X,)=s0,  where  X  is  the  average  number  of  event  per  time  interval.  As  stated,  our  process  of 
seismic  event  occurrence  is  non-Poissonian  but  we  hope  that  the  same  condition  should  be  met  to  use  the  asymptotic 
distribution  for  maxima.  30  events  per  interval  gives  the  value  exp(-?i)  <  10'^"^ . 
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Figure  5.3  Probability  distribution  of  the  20  days  maximum  seismic  energy  release  at  the 
coal  level  501 
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Figure  5.4  Probability  distribution  of  the  20  days  maximum  seismic  energy  release  at  the 
coal  level  504 
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Figure  5.5  Probability  distribution  of  the  20  days  maximum  seismic  energy  release  at  the 
coal  level  507 
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Figure  5.6  Probability  distribution  of  the  20  days  maximum  seismic  energy  release  at  the 
coal  level  510 
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complex  structure  with  two  more  high  energy  modes.  The  limiting  value  of  the  first  from  these 
modes  is  located  in  the  second  half  of  order  10^  J,  and  that  of  the  second  one  exceeds  10*  J. 

The  distributions  of  energy  maxima,  associated  with  the  coal  levels,  were  compared  by 
means  of  the  Kolmogorov-Smimov  two  sample  test.  The  results  of  the  test  are  given  in  Table 
5.1 


Table  5.1.  Results  of  Kolmogorov-Smirnov  test  for  distributions  of  maximum  seismic  energy  released  at  tbe  coal 
levels 


Coal  level 

Naof 

observations 

Significance  level 

416 

501 

504 

507 

510 

416 

28 

X 

<0.001 

>0.05 

>0.05 

<0.01 

501 

162 

X 

<0.001 

0.001 

<0.001 

504 

14 

X 

>0.05 

>0.05 

507 

21 

X 

>0.05 

510 

50 

X 

The  test  does  not  point  out  significant  differences  in  the  groups  comprising  the  samples 
associated  with  the  levels:  416,  504,  507  and  comprising  the  samples  associated  with  the  levels: 
504,  507,  510.  The  level  501  associated  sample  tumes  out  to  be  highly  different  from  all  other 
samples  connected  with  all  other  levels.  The  similarities  in  the  sanqtle  distributions  form  groups 
of  vertically  ordered  levels.  Thus  it  is  likely  that  the  process  for  enei^  maxima  changes 
gradually  with  increasing  depth  of  exploitation.  The  level  501,  outstanding  from  the  point  of 
view  of  general  rockbursts  danger,  tumes  out  to  be  also  outstanding  regarding  the  maximum 
energy  generation. 

The  probability  distributions  of  maximum  seismic  energy  release  at  the  particular  minp. 
region,  in  the  studied  period,  are  shown  in  Figures  5.7  -  5.11,  respectively.  The  following 
modes  of  the  distributions  were  identified: 

-NWpart:  I  -  E^5xl0^  J  ,  H  -  E^«10*  J  ,  UI  -  E„^«5xl0' J  ,  IV-E^«3xl0"j 
and  possibly  V  -  En^>  10^  J  witii  (Mily  one  observation; 

-SWpart:  possibly  I  -  E^«5xl0^  J  ,  II-E^«6xlO^J  ,  m-E^«10^J  , 

IV-  E^xlO'j  ,  V-E^>5xl0’j; 
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EMPIRICAL  DISTRIBUTION  OF  ENERGY  MAXIMA 
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-Cpart:  I-E^lO'j  ,  n-E^«10"j  ,  m-E^«10’j; 

-  NE  part:  I  -  En^ax^SxlO'’  J  ,  II  -  E^  >  10^  J  with  only  one  observation; 

-SEpart:  I-E^«5xl0'j  ,  n-E^^lO’j  ,  III  -  Ej„^»2xl0*  J  ,  IV-Ew“10’j. 
Like  previously,  these  distributions  of  maxima  were  compared  by  means  of  the 
Kolmogorov-Smimov  two  sample  test.  The  results  are  given  in  Table  5.2. 


Table  5.2.  Results  of  Kolmogorov-Smirnov  test  for  distributions  of  maximum  seismic  energy  released  at  the 


mine  regions 

Mine  region 

No.  of 

observations 

Significance  level 

NW  SW  C 

NE 

SE 

INW 

28 

X  <0.001  <0.01 

>0,05 

<0.001 

sw 

162 

X  <0.001 

<0.001 

<0.05 

504 

14 

X 

<0.05 

>0.05 

507 

21 

X 

<0.01 

510 

50 

X 

The  distributions  of  logE  maxima  of  the  seismic  series  associated  with  NW  and  NE  parts  of  the 
mine  and  the  distributions  of  logE  maxima  associated  with  C  and  SE  parts  are  similar  (not 
significantly  different). 

The  picture  provided  by  this  study  of  the  regions  associated  distributions  of  logE  maxima 
is  less  clear  than  that  obtained  from  the  previous,  level  associated  samples  study  .  The  number  of 
inflexions  in  the  distributions  increased,  hence  the  estimates  of  the  limiting  values  are  less 
reliable.  Doubtlessly  the  studies  prove  that  there  are  different  processes  of  event  generation.  It 
seems  that  they  have  different  dynamics  in  different  parts  of  the  mine  because  the  limiting  values 
of  the  distributions  of  maxima  are  different.  We  can  join  these  processes  in  crudely  determined 
three  classes.  Class  one  of  low  energy  mode  jM'ocesses  contains  modes  I  and  H,  with  limiting 
values  in  the  second  half  of  the  order  10“^  J  or  m  the  first  half  of  the  order  10*  J,  These  are 
probably  processes  directly  controlled  by  mining  works  in  particular  excavations,  taking  place  in 
the  coal  layers  or  in  their  direct  surrounding.  The  medium  energy  modes,  with  limiting  values 
about  10^  J  assigne  the  second  class.  The  second  class  processes  are  probably  taking  place  in 
thick  sandstone  layers  separating  the  coal  seams  and  are  controlled  by  the  stress  field  resulted 
integrated  effects  of  the  active  excavation,  old  workings  geometry  and  local  tectonic  features. 
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The  third  class,  which  comprises  high  energy  mode  processes  could  be  linked  to  generalized 
regional  stress  field  variations. 

We  are  unable  to  give  the  detailed  explanation  of  the  processes  of  seismic  event 
generation  upon  the  extreme  value  statistics  study  alone.  Certainly  the  generating  mechanisms 
are  more  complex  than  these  mentioned  here.  Their  better  recognition  needs  complementing  the 
statistical  analysis  with  the  rock  mechanics  modeling  and  the  physics  of  rockmass  ftacturing 
studies. 

Conclusions 

1.  The  maximum  energy  value  statistics  method  tumes  out  to  be  a  valuable  tool  in  studying 
complexity  of  the  series  of  seismic  events. 

2.  The  results  of  the  study  prove  that  the  energy  distribution  of  seismic  events  induced  by 
nuning  is  multimodal  in  both  mine-wide  scale  as  well  as  at  particular  coal  layers  and  in 
particular  mine  parts. 

3.  The  maximum  value  studies  enable  to  single  out  particular  modes  of  event  energy 
distributions  and  to  estimate  roughly  their  limiting  values.  It  is  expectable  that  different  event 
generating  processes  are  responsible  for  the  different  modes  of  energy  distributions. 

4.  Some  coal  levels  and  some  mine  parts  have  similar  distributions  of  maximum  released 
energy,  other  ones  have  significantiy  different. 

5.  The  distributions  of  energy  of  events  recorded  in  strongly  seismically  active  mine  parts  or  at 
strongly  active  coal  levels  have  more  high  energy  modes  than  the  distributions  associated  with 
mine  regions  where  the  seismic  activity  was  moderate. 

6.  In  general  three  classes  of  processes  of  seismic  event  generation  are  identified.  The  processes 
belonging  to  these  classes  are  responsible  for  low,  medium  and  high  energy  modes  of  event 
energy  distribution,  respectively.  The  limiting  values  of  these  processes  depend  upon  the  part  of 
the  mine  they  take  place. 

7.  Determining  of  controlling  factors  and  mechanisms  of  the  identified  processes  will  not  be 
possible  unless  information  firom  the  rock  mechanics  and  physics  of  fracturing  studies 
complements  results  of  the  statistical  analysis. 
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8.  The  extreme  value  statistics  is  an  indirect  way  to  identify  processes  of  seismic  event 
generation.  It  is  not  possible,  on  a  basis  of  this  method,  to  separate  from  the  seismicity  series 
events  generated  by  specific  identified  processes. 
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6. 


Deflections.  A  preliminary  analysis. 


It  is  well  known  that  series  of  mining  induced  events  are  not  distributed  randomly  in 
space  but  have  a  certain  geometrical  time-varying  structure.  The  seismic  events  tend  to  occur 
around  plane  structures  like  old  workings  edges,  faults  and  other  weak  zones.  Various  methods 
have  been  applied  to  analyze  geometrical  distribution  of  sources  in  order  to  group  events  that  are 
due  to  a  single  ftacturing  process.  These  are  cluster  analysis  techniques  where  relations  between 
events  are  inferred  from  their  mutual  distances,  in  a  hope,  that  different  generating  processes 
would  result  in  creating  different  clusters  or  different  groups  of  clusters  of  induced  seismicity 
(e.g.,  Frohhch  and  Davis  1986,  Davis  and  Frohlich  1991,  Kijko  et  al.  1993).  According  to  our 
experience  gained  from  trials  made  on  the  seismic  catalog  from  Wujek  mine,  the  clustering 
methods  may  work  well  when  the  catalog  is  of  small  or  moderate  size.  For  large  series,  and 
events  which  are  described  by  only  four  -  five  parameters,  the  cluster  classification  leads  to 
either  trivial  or  non-unique  results,  difficult  to  interpret.  It  is  also  not  still  resolved  how  to  define 
properly  the  distance  in  the  time-space,  to  enable  four  dimension  clustering,  not  mention  about  a 
higher  order  space  like  the  time-location-energy  one.  We  evidenced  formerly  that  the  mining- 
induced  seismicity  is  nonstationary  in  both  time  and  space.  Hence  significant  clustering  is 
expected  rather  in  four,  five-dimensional  than  in  only  the  Euchdean  space. 

Geometrical  approach  is  also  represented  by  the  fractal  analysis  of  seismic  source 
distribution  (e.g.,  Smalley  et  al.  1987,  Xie  and  Pariseau  1992,  Idziak  and  Zuberek  1995, 
Lasocki  and  Mortimer  1996).  Variations  of  the  fractal  dimensions  among  the  seismic  series, 
when  the  series  have  fractal  character,  or  values  of  different  order  fractal  dimensions,  in  case  of 
multifiactal  structure  of  the  series,  are  tried  to  be  interpreted  in  terms  of  various  generating 
processes.  Still  studies  are  needed  concerning  properties  of  the  fractal  dimension  distributions 
which  could  enable  distinguishing  of  the  fractal  dimensions  (Cosentino  et  al.  1996).  The  other 
drawback  of  the  fractal  methods  is  the  same  as  connected  with  the  cluster  analysis:  the  problem 
with  defining  distance  in  the  multidimensional  space. 

Limited  number  of  parameters  of  events  limits,  in  our  opinion,  possibilities  to  study  the 
complex  structure  of  the  seismic  series.  Therefore  we  have  introduced  and  investigated  a  new 
parameter  to  describe  the  directional  structure  of  the  series.  The  trend  of  the  seismic  series  is 
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characterized,  in  our  study,  by  the  distribution  of  deflection  of  straight  lines  connecting  the 
epicenters  of  every  two  consecutive  events.  The  deflection  is  measured  from  the  NS  direction. 
The  sense  of  the  vector  connecting  epicenters  is  not  assessed  and  the  identity 


a +  180° 
a -180° 


for  a  e[-180°,-90°] 
for  a  e  [90°  ,180°] 


(6.1) 


is  set.  In  this  way,  the  deflection  varies  in  the  range  of  [-90°  ,  90°  ]. 

The  deflection  is  regarded  as  a  random  variable  with  a  certain  probability  distribution. 
Modes  of  the  probability  distribution  function  of  deflections  show  dominant  trends  of  the  studied 
seismic  series.  The  parameter  forms  a  time  series,  hence,  in  addition  to  information  about 
geometrical  assembly  of  epicenters,  it  kept  also  track  of  time  ordering  of  events.  Since,  as  it  was 
mentioned,  the  plane  structures  (linear  on  the  horizontal  plane)  are  likely  to  play  the  major  role 
in  generating  seismic  events,  we  assume  that  the  events  building  the  specific  dominant  trend  are 
outcomes  of  a  single  generating  process. 

For  the  N-event  series  we  receive  the  (N-l)-point  set  of  deflections.  Thus,  unlike  in  the 
maximum  value  studies,  the  analysis  of  deflection  can  be  performed  on  particular  seismicity 
clusters.  The  clusters  of  seismicity  at  Wujek  mine  were  identified  in  Section  3.  The  descriptive 
statistics  studies  showed  that  these  clusters  differ  in  statistical  properties.  Therefore,  the  analysis 
of  deflections  was  made  on  the  selected  clusters.  In  cases  when  the  sample  size  was  large  enough 
for  such  the  analysis,  all  complete  series  as  well  as  their  subsets,  built  by  removing  events  below 
the  prescribed  energy  thresholds,  were  considered.  The  studied  data  sets  are  presented  in  Table 


6.1 


The  sets  1-3,  5-12,  14-19,  21-23,  16-30  were  built  form  the  seismic  series  that  occurred 
in  vicinity  of  active  excavations  (see:  Figs.  3.3-3.17).  Other  sets  concern  cases  of  "strange” 
activity,  which  caimot  be  linked  to  mining  works.  C  part,  at  the  level  501,  was  divided  with  the 
W-S  line  at  x=  18700  into  the  southern  area  with  mining  stopes  and  the  northern  area  where  no 
exploitation  could  justily  occurrence  of  seismic  events  (see:  Fig.  3.6).  Also  NE  part  of  the  level 
510  was  split  with  the  N-S  line  at  y=-9850,  and  only  its  western  part,  that  had  no  obvious 
connection  with  mining,  was  analyzed  (see:  Fig.  3.15). 
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Table  6. 1  Data  sets  studied  in  the  analysis  of  deflections 


No' 

Coal  level 

Mine  part 

Energy 

threshold 

Number  of 

events 

1 

416 

C 

no 

1342 

2 

416 

C 

3.2xl0’j 

738 

3 

416 

c 

lO"! 

134 

4 

416 

sw 

no 

89 

5 

501 

sw 

no 

4816 

6 

501 

sw 

lO^J 

1207 

7 

501 

sw 

3.2xl0^J 

344 

8 

501 

sw 

icfj 

49 

9 

501 

southern  part  of  C 

no 

1572 

10 

501 

southern  part  of  C 

lO^J 

524 

11 

501 

southern  part  of  C 

3.2xlO^J 

158 

12 

501 

southern  part  of  C 

10*J 

49 

13 

501 

northern  part  of  C 

no 

112 

14 

504 

C 

no 

367 

15 

504 

C 

3.2xl0’j 

197 

16 

504 

c 

lO^J 

71 

17 

507 

NW 

no 

1199 

18 

507 

NW 

3.2xl0’j 

889 

19 

507 

NW 

lO^J 

314 

20 

507 

SW 

no 

49 

21 

510 

NW 

no 

272 

22 

510 

NW 

3.2xl0’j 

168 

23 

510 

NW 

lO^J 

52 

24 

510 

western  part  of  NE 

no 

357 

25 

510 

western  part  of  NE 

3.2xl0’j 

230 

26 

510 

western  part  of  NE 

10"J 

56 

27 

510 

SE 

no 

1677 

28 

510 

SE 

3.2xlO^J 

1133 

29 

510 

SE 

lO^J 

348 

30 

510 

SE 

3.2xlO^J 

91 

For  each  data  set  the  respective  set  of  deflections  was  prepared.  Inspection  of  the 
histograms  of  deflections  for  the  considered  seismic  series  suggests  their  nonuniform 
distribution.  At  least  one  distinct  mode  of  the  distribution  was  visible  in  each  histogram.  Many 
histograms  had  also  the  secondary  modes.  In  some  cases  neither  the  main  nor  secondary  modes 
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were  located  at  the  same  places  for  different  subsets  of  the  same  data  sets.  Figure  6. 1  shows,  as 
an  example,  the  histograms  of  deflections  for  sets  3  and  5,  that  is  for  the  same  seismicity  cluster 
from  SW  part  of  the  level  501,  with  no  energy  threshold  and  the  threshold  of  lO^J  applied, 
respectively.  The  main  mode  of  the  distribution,  located  at  some  50°  for  set  3,  is  shifted  to  some 
65°  for  set  5.  The  strong  event  series  (set  5)  shows  also  the  second  trend  of  some  -20°. 

In  order  to  verily  the  suggestion  about  the  nonuniform  distribution  of  deflections,  the 
empirical  distributions  were  compared  with  the  theoretical  uniform  distribution  by  the  chi-square 
and  one-sample  Kolmogorov-Smimov  tests.  For  all  subsets  both  testing  procedures  have 
rejected,  with  high  significance,  the  hypothesis  about  uniform  character  of  the  empirical 
distributions  of  deflections. 

In  practically  all  analyzed  cases  the  histograms  of  deflections  imply  that  the  distribution 
of  deflections  of  seismic  events  has  complex  multimodal  structure.  The  identification  of  the 
modes  is  supposed  to  deliver  valuable  information  about  components  of  the  data  and,  indirectly, 
about  processes  involved  in  the  seismicity  generation.  Because  of  the  subjective  choice  of  a  way 
of  binding  the  histogram  of  random  sample  is  always  a  subjective  structure.  The  more  complex 
and  multimodal  distribution  of  the  variable  is,  the  smaller  sample  is  analyzed,  the  less  certain 
shape  of  the  distribution  is  assessed  from  the  histogram.  Although,  in  general,  our  study 
concerns  large  size  samples  we  intend  to  work  also  with  sets  of  some  50  or  less  points  where  a 
final  shape  of  the  histogram  will  dominantly  depend  upon  the  way  of  binding.  We  need, 
therefore,  methods  which  can,  more  objectively  than  the  histograms  do,  recover  probability 
density  functions  from  random  samples.  The  nonparametric  kernel  estimation  of  probability 
density  function,  the  technique,  so  far  as  we  know,  for  the  first  time  used  in  studies  of 
seismicity,  seems  to  fit  to  our  purposes.  A  theoretical  development  of  the  method  and  results  of 
various  experiments  necessary  before  applying  the  method  to  real  data  are  given  in  the  next 
section. 


NO  OF  OBS  NO  OF  OBS 


Summary  and  conclusions 

1.  A  new  way  to  study  nonrandom  character  of  the  spatial  distribution  of  seismic  series  was 
proposed.  The  method  is  based  on  the  deflection  of  straight  lines  connecting  epicenters  of  eveiy 
two  consecutive  events,  measured  from  the  NS  direction. 

2.  Distributions  of  deflection  for  seismic  series  from  Wujek  coal  mine  turned  out  to  be  complex 
and  not  uniform.  There  were  clear  dominant  and  sometimes  also  secondary  directions.  The 
modes  of  the  distributions  are  expected  to  be  formed  by  events  originated  by  different  generating 
processes. 

3.  The  identification  of  modes  of  deflection  distributions  cannot  be  made  from  histograms 
because  the  histograms  strongly  depend  on  a  subjective  choice  of  binding.  As  an  alternative 
approach  the  nonparametric  kernel  estimation  of  the  probability  density  function  is  proposed. 
We  think  that  this  method  can  be  used  not  only  in  the  analysis  of  deflection  but  also  to  study 
structure  of  the  seismic  energy  or  logarithm  of  energy  distributions. 
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7.  Nonparametric  kernel  estimation  of  probability  density  functions 


Statistical  inference  on  a  probability  density  function  (p.d.f.),  given  random  sample,  is 
commonly  made  by  assuming  a  function,  estimating  its  parameters  and  assessing  its  validity 
using  an  adequate  test.  Then,  for  visual  inspection  rather  than  for  other  purposes,  the 
estimated  probability  density  function  can  be  compared  to  a  histogram.  Unlike  empirical 
probability  distribution  function  or  empirical  probability  exceedance  function,  histogram  is 
a  highly  subjective  strucmre  and  its  resemblance  to  the  true  p.d.f.  can  always  be  criticized. 
Nevertheless,  there  is  permanent  need  for  developing  methods  enabling  the  user  to  construct 
a  reliable  estimator  of  a  p.d.f. ,  especially  where  rationale  exists  that  the  true  p.d.f.  is  a  multi¬ 
peak  function.  Of  course,  a  p.d.f.  can  be  derived  from  a  given  cumulative  probability 
distribution  function  (cdf).  Unfortunately,  the  latter  is  an  integral  function,  which  means  that 
it  screens  local  effects  and  insufficiently  large  peaks  are  not  detectable. 

Nonparametric  kernel  estimation,  started  with  papers  of  Rosenblatt  (1956)  and  Parzen 
(1962),  is  an  approach  which  deals  with  direct  estimation  of  unknown  true  p.d.f.  using 
random  sample  data.  Up  to  now,  many  estimators  have  been  developed,  their  theoretical 
properties  recognized.  However,  from  practical  point  of  view,  the  most  important  question 
has  not  been  sufficiently  answered  yet:  How,  based  on  a  random  sample  of  real  data,  is  to 
be  selected  the  smoothing  parameter  hi.  Because  theory  cannot  and  will  not  fully  answer  to 
this  question,  the  only  reasonable  way  is  the  experiment.  It  implies  two  possibilities 
complementing  each  other:  (a)  the  Monte  Carlo  experiments  enabling  the  researcher  to 
investigate  many  features  which  he  wants  to  know,  and  (b)  applications  of  nonparametric 
methods  to  real  data  samples  and  using  the  physical  information  underlying  the  samples  to 
verify  (i.e.  to  confirm  or  reject)  the  statistical  findings  discovered  by  the  nonparametric 
method  used.  However,  although  the  possibility  (a)  is  rather  widely  employed,  real  data  are 
used  for  illustration  purposes  rather  than  as  an  additional  source  of  information  which 
qualitatively  differs  from  generated  pseudorandom  numbers  and,  for  this  reason,  is  the 
information  of  special  importance. 

Our  work  is,  on  the  one  hand,  connected  with  the  possibility  (a),  when  we  illustrate  how 
kernel  methods  work  in  cases  which  are  of  our  interest,  and,  on  the  other  hand,  it  also 
employs  the  possibility  (b)  in  order  to  identify  the  real  distributions  of  deflections  and 
logarithms  of  energy  of  tremors,  as  well  as  in  order  to  show  how  the  relevant  physical 
mformation  correlates  with  the  approximations  obtained  and  justifies  their  use.  It  should  be 
stressed  here  that  while  numerous  papers  have  been  published  on  theoretical  and/or  Monte 
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Carlo  findings,  papers  dealing  strictly  with  application  problems  of  the  p.d.f.  nonparametric 
kernel  estimation  are  rare  and  to  our  best  knowledge  this  approach  has  not  been  yet  used  in 
seismicity. 


Kernel  estimation  of  probability  density  function 


A  kernel  estimator  j^x)  of  a  p.d.f.  fix)  can  be  expressed  (Parzen  1962)  as 


hx)  =  -L'ts 

nhU 


x-x. 

~ir 


(7.1) 


where  x„  i = 1 ,2, . . .  ,n,  is  a  random  sample,  K(.)  is  a  given  kernel  function  and  A  is  a  positive 
smoothing  parameter  (the  terms  bandwidth  or  window  width  are  also  used).  Under  conditions 
given  in  Parzen’s  (1962)  paper,  the  j^x)  is  asymptotically  unbiased  and  consistent.  For  the 
remainder  of  this  section  it  is  assumed  that  the  sample  x^,  /=l,2,...,n,  is  nondecreasing. 

In  practice,  both  a  kernel  function  and  smoothing  parameter  have  to  be  chosen.  The  choice 
of  a  kernel  function  is  not  very  important  and  many  symmetric  unimodal  functions  can  give 
similar  efficiency  of  the  method.  For  example,  Silverman  (1986)  compared  the  efficiency 
eff(K)  of  five  following  symmetrical  kernel  functions  K(y): 
a  parabola: 


K(y)  =  . 

3 

fi-z!! 

for  |y|  <  5^^ 

4v^ 

5 

eMK)  =  1 

0 

for  |y|  >5^'^ 

a  parabola  of  fourth  order; 


(7.2) 


a  triangle; 


K(y) 


for  |y|  <  1 
0  for  \y\>l 


ejfiK)  =  0.9939 


K(y) 


Ji-bl 

1  0 


for 

for 


y  <1 

y  >1 


eMK)  =  0.9859 


(7.3) 


(7.4) 


a  rectangle: 


^0-) 


1/2 

0 


for 

for 


y  <1 

y  >1 


eMK)  =  0.9295  (7-5) 
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and  a  Gauss: 


=  —  exp 
sjlTT 


ejfiK)  =  0.9512  (7.6) 


Thus  practically  each  symmetric  unimodal  function  can  be  applied  as  a  kernel  function. 

The  choice  of  h  is  more  crucial  for  performance  of  the  method  and  several  procedures 
exist  for  estimation  of  the  smoothing  parameter  value  (eg.:  Scott  and  Factor  1981,  Katkovnik 
1985,  Scott  and  Terrel  1987):  method  of  minimization  of  integral  mean  square  error,  cross- 
validation  (modified  maximum  likelihood  function)  method,  k-nearest  neighbor  method,  a 
method  based  on  the  Fourier  series  expansion  of  a  true  density  (in  general:  orthogonal  series 
estimators),  and  others.  According  to  Scott  and  Factor  (1981)  none  of  them  seems  to  be  better 
than  the  others  so  the  authors  suggest  constructing  several  estimators  in  order  to  eliminate 
occasional  poor  ones. 


Accuracy  measures  of  p.d.f.  estimators 

The  accuracy  of  the  kernel  estimator  (7.1)  can  be  expressed  in  many  ways  of  which  the 
mean  square  error,  MSE,  and  integral  mean  square  error,  IMSE,  are  applied. 

The  mean  square  error,  MSE,  is  an  error  estimation  at  x: 

MSE{x,h)  =  E\f{x)-f{x)f  (7.7) 

It  can  be  decomposed  into  two  basic  components: 

MSE  =  a^(f)  ^  bias^(^  (7-8) 

where  (^(f)  is  a  random  error  component: 

<7^(/)  =  E{f(x)-Ef(x)]^  (7.9) 

and  bias^(f)  is  a  systematic  error  component: 

bias^(J)  =  E[Ef(x)-f(x)f  =  \Ef(x)-f(x)\^  (7-10) 

If  /z  -»  0  then  l^x)  -»■  f(x)  so  the  bias  tends  to  zero.  However,  the  random  error  o(^ 
increases  then  to  infinity  because  the  estimator  /tends  to  a  linear  combination  of  Dirac  delta- 
functions.  Thus  the  choice  of  the  smoothing  parameter  h  (an  optimum  value  of  h)  should  be 
a  compromise  between  these  two  errors. 
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The  integral  mean  square  error,  IMSE,  is  an  estimation  of  the  global  error  for/: 

IMSE  =  j  E\S{x)-J{x)]'^dx  (7.11) 

Minimization  of  the  IMSE  as  a  function  of  h  is  one  of  the  most  commonly  used  techniques 
for  estimation  of  the  smoothing  parameter  h. 

Methods  for  estimation  of  smoothing  parameter 


The  choice  of  h  is  the  most  crucial  for  the  performance  of  any  method.  If  h  is  too  small, 
the  resulting  estimated  density  function /(.x)  is  too  little  smoothed  and  shows  too  much  noise 
contained  in  the  sample.  And  on  the  contrary:  if  h  is  too  large,  the  resulting  oversmoothing 
masks  the  structure  of  the  true  p.d.f.  fiX).  Unfortunately,  the  choice  of  each  method  is  to 
some  extent  subjective,  and  it  is  upon  a  researcher  to  decide  which  method  is  to  be  applied. 
The  situation  seems  to  be  similar  to  the  problem  of  estimating  the  spectral  density  function, 
especially  if  the  data  is  not  very  numerous. 

We  have  selected  and  decided  to  compare  two  of  the  existing  methods.  The  first  is  a 
simplification  of  the  IMSE  minimization  method,  and  the  second  one  is  based  on  the  concept 
of  modified  (empirical)  likelihood  function.  Anticipating  a  little,  it  can  be  said  that  the  first 
method  is  the  simplest  and  most  reliable  because  is  simple  for  calculation  and  gives  one  (and 
only  one)  h  value  for  every  sample.  The  second  method  is  more  time  consuming,  usually 
gives  more  than  one  h  value,  and,  not  infrequently,  h=0. 


The  optimum  value  of  the  smoothing  parameter  h  can  be  obtained  by  minimization  of  the 
IMSE.  According  to  Scott  and  Factor  (1981),  this  leads  to  the  following  asymptotical  optimal 
choice  for  h  : 


h  =  ci{K)  J3(f)  n-^'^ 


(7.12) 


where  ci{K)  and  /3(/)  are  factors  depending  on  kernel  function  K  and  p.d.f.  /  respectively: 


11/5 

oi(K)  =  ^Kiyfdy 


-2/5 


I  K(y)y^dy 


m 


(7.13) 

(7.14) 


The  unknown  function /in  the  expression  /?(/)  is  usually  replaced  with  its  estimator  based  on 
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the  chosen  kernel  function  K.  The  method  has  good  properties,  as  it  was  numerically 
demonstrated  that  the  choice  of  h  based  on  (7.13)  is  in  the  average  nearly  optimal  for  even 
as  small  samples  as  n=25  (Scott  and  Factor  1981). 

It  is  possibly  to  obtain  simpler  results.  It  was  shown  (Adamowski  and  Feluch  1987,  Feluch 
1994)  that  for  the  parabolic  kernel 


K(y) 


3 

4^/5^ 


0 


for  |y|  <^f^ 
for  |y|  >V^ 


(7.15) 


the  smoothing  parameter  h  can  be  approximated  by  an  expression 


h  < 


5^(2/-«-1)x,. 

/  =  1 


n 


(7.16) 


Because  the  value  of  b  does  not  influence  the  error  of  the  estimator  (7.1),  the  authors 
propose  to  use  the  value  b=ll5,  which  leads  to  simpler  equations  for  the  kernel  K 


K(y) 

and  the  smoothing  parameter  h 


l(l-y^)  for  |y|  <  1 
0  for  \y\>l 


J^(2i-n-l)x,. 

1  =  1 


(7.17) 


(7.18) 


The  sign  of  the  weak  inequality  in  (7. 16)  has  been  replaced  in  equation  (7. 18)  by  the  sign  of 
equality.  The  h  value  thus  obtained  is  the  smallest  upper  limit  for  smoothing  parameter 
produced  by  the  method  (7.12)  of  minimizing  IMSE  (7.11). 


Another  method  of  obtaining  a  smoothing  parameter  estimator  is  based  on  the  modified 
(empirical)  likelihood  function  L{h)  (Scott  and  Factor  1981,  Katkovnik  1985,  Silverman 
1986): 
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m  =  n/  (•^/) 


(7.19) 


in  each  case  and  is  simple  and  fast. where  is  equal  to 


=  -1  K  ^ 


(7.20) 


f^\xj)  is  an  estimator  of  the  p.d.f.  at  Xj  based  on  the  whole  sample  excluding  point  Xj. 
Combining  the  last  two  equations,  we  get 


L{h)  =  n 


(n-l)/z,=tT 


(7.21) 


The  optimum  value  h  of  his  computed  from  the  condition: 


h  =  argmaxL{h) 

h&O 


(7.22) 


The  procedure  is  quite  effective  even  for  small  n  (Katkovnik  1985).  However,  when  kernel 
(7.17)  is  used,  possible  values  of  h  must  be  greater  than  /z^in 


^min  =  maxj  Ax^,  mm(Ax,.,Ax,..i),  Ax  I,  Ax,.  =  x.-x,.., 


(7.23) 


otherwise  L{h)  =  0.  The  method  is  thus  insensitive  to  h  values  from  the  interval  (0,h^i„). 


Monte  Carlo  visualization  of  the  f(x)  estimator  aualir 


Unlike  parametric  methods  of  distribution  function  estimation,  nonparametric  methods 
employ  local  estimation,  not  global  one.  The  resulting  estimator  of  p.d.f.,  ^x),  is  basically 
assigned  to  that  interval  Ax  of  realizations  of  random  variable  X  which  includes  the  random 
sample  Xj,  i=l,2,...,n,  Ax=:x„-Xi,  not  covering  thus  the  whole  interval  of  realizations  ofX 
Extrapolation  towards  low  or  high  values  of  X  is  possible  for  only  a  limited  left  vicinity  of 
Xj  and  limited  right  vicinity  of  x„;  as  far  as  the  kernel  and  smoothing  parameter  make  it 
possible.  Methods  for  kernel  estimation  quality  assessing  are  also  different  from  those  for 
parametric  estimation  and  consist  mainly  of  assessment  of  the  theoretical  quality  of  the 
smoothing  parameter  estimator  (consistency,  bias,  etc)  and  the  quality  assessed  from 
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simulation  experiment,  a  Monte  Carlo  simulation.  It  enables  the  researcher  to  observe  the 
interesting  aspects  of  the  studied  problem  and  find  its  pros  and  contras.  The  Monte  Carlo 
simulation  results,  which  are  presented  on  subsequent  pages,  are  to  illustrate  some  main 
features  of  the  kernel  estimation  method  for  p.d.f.  estimation.  Throughout  this  section,  the 
smoothing  parameter  h  is  calculated  by  equation  (7.18)  which  is  the  simplest  method,  giving 
somewhat  overestimated  value  (with  the  resulting  oversmoothing). 

In  order  to  show  various  aspects  of  the  kernel  estimation  of  a  p.d.f.,  the  Monte  Carlo 
simulation  was  applied  to  some  unimodal  and  bimodal  distributions.  We  selected  those 
distributions  that  resembled  the  shape  of  empirical  event  energy  or  logarithm  of  energy 
distributions  or  possessed  features  expected  to  be  met  in  distribution  of  deflection.  In  each 
case  a  true  p.d.f.  is  assumed  to  be  known,  a  random  sample  was  taken,  and  the  kernel 
estimation  was  used  to  produce  an  estimator  of  f(x). 


Lognormal  unimodal  distribution 


This  two-parameter  distribution,  defined  by  the  following  equation 


yi„(^;M>cf) 


1 

- —exp 

OXyJl'K 


(lnx-/i)M 


2a^ 


(7.24) 


was  used  in  its  standard  version  (/^=0  and  <7=1)  to  produce  18  random  samples  of  size  n 
equal  to  30,  60,  90,  120,  150  and  180,  each  size  case  sampled  three  times.  Results  are 
presented  in  Fig.  7.1  where  also  the  theoretical  p.d.f.  is  drawn.  The  maximum  for  that 
function  is  reached  at  the  point  x  =  expiix-a^)  =  e  «  2.7. 

Lomnitz  (1974)  postulated  the  use  of  lognormal  distribution  to  model  an  empirical 
distribution  of  earthquake  magnitudes.  The  distribution  was  also  tried  to  represent  logarithm 
of  seismic  energy  in  mining  induced  seismicity  (Lasocki  1989). 

Generally,  when  looking  down  the  Fig.  7.1  along  with  the  increasing  sample  size  n,  the 
expected  tendency  is  visible.  The  quality  (resemblance  to  the  theoretical  p.d.f.)  of  simulated 
p.d.f.-s  seems  to  increase  with  increasing  n.  In  certain  cases,  for  small  sample  sizes,  there 
exist  some  qualitative  differences  as,  for  example,  distinctively  developed  bimodality  or  shifts 
of  estimated  maxima  of  J{x)  when  compared  to  the  true  p.d.f.  f(x).  As  n  increases,  the 
differences  tend  to  be  smaller:  the  bimodality  is  disappearing,  locations  of  the  true  and 
estimated  maxima  are  closer  one  to  another,  and  the  shapes  of  estimated  f(x)  and  true 
p.d.f.  f(x)  are  more  consistent,  independently  to  the  sample  range  (e.g.  0  -  60  for  the  sample 
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Fig.  7.1.  Lognormal  distribution  p,d.f.-s:  theoretical  (equation  (7.24)),  ^—0,  a— 1,  (heavy  line)  and  kernel  estimated 
(thin  line  with  circles)  for  various  samples  of  the  same  size  n  (rows)  and  various  sizes  (columns);  h  in  the  n/h 
expression  is  the  smoothing  parameter  calculated  by  equation  (7.18). 


180/0.63  or  0  -  30  for  120/0.55  in  Fig.  7.1).  In  each  case  the  estimated  function  resembled 
the  true  p.d.f.  at  its  both  tails  well.  This  is  the  result  of  the  continuous  falling  to  zero  of  both 
tails  of  the  theoretical  p.d.f.  It  would  be  interesting  to  see  how  the  discontinuity  of  the  true 
p.d.f.  affects  its  estimator  properties.  For  this  purpose  the  truncated  exponential  distribution 
properties  are  investigated  in  the  next  subsection. 

Truncated  exponential  distribution 

The  one-parameter  truncated  exponential  distribution  is  defined  in  the  (0,1)  interval  as 
follows: 


f,eM 


l-e-“ 


0<x<  1 


0  otherwise 


(7.25) 


Corrections  introduced  to  the  Gutenberg-Richter’s  relation,  that  accounted  for  a  deficit  of 
the  number  of  events  in  the  large  magnitude  range,  led  to  the  truncated  exponential 
distribution  for  earthquake  magnitudes  (Cosentino  et  al.  1977).  At  present,  the  truncated 
exponential  distribution  is  the  most  widely  used  model  for  logarithms  of  energy  in  strong 
tremor  hazard  analyses  of  induced  dynamic  failure  phenomena  (Lasocki  1992a,  1993, 
Gibowicz  and  Kijko  1994,  Kijko  and  Funk  1994). 

The  p.d.f.  is  not  continuous  at  two  points:  x=0  and  x=\.  This  discontinuity  causes 
additional  problems  in  reproducing  the  true  p.d.f.  f{x)  as  the  edge  effect  occurs.  Its  influence 
on  the  estimated  p.d.f.  can  be  strong,  especially  when  the  jumps  at  the  discontinuity  points 
are  large,  and  can  be  seen  in  Fig.  7.2  as  the  tendency  towards  smooth  transition  between  zero 
and  non-zero  values  of  an  estimated  p.d.f.  The  graphs  in  Fig.  7.2  shows  that  this  tendency 
is  present  in  all  sampled  distributions  and  it  seems  that  the  increase  from  n=30  to  n=180 
does  not  improve  the  resemblance  of  estimates  ^x)  to  the  true  f{x)  at  the  right  vicinity  of 
point  x=0  and  the  left  vicinity  of  point  x= 1.  Multimodality  is  stronger  than  for  the  previous 
(lognormal)  distribution  because  the  edge  effect  supports  strongly  the  production  of  additional 
modes  of  the  estimated  p.d.f.  This  can  be  seen  especially  around  the  left  border  where  the 
p.d.f.  values  are  about  three  times  greater  then  for  the  right  one.  Central  parts  of  the 
estimated  p.d.f.-s  seem  to  be  less  affected  by  the  edge  effect,  which  is  visible  in  Fig.  7.2  as 
decreasing  amplitudes  of  estimated  p.d.f.  fluctuations  around  the  true  p.d.f.  as  the  sample  size 
n  increases. 
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Fig.  7.2.  Truncated  exponential  density  functions:  theoretical  (equation  (7.25)),  a  =  l,  (heavy  line)  and  kernel 
estimated  (thin  line  with  circles)  for  various  samples  of  the  same  size  n  (rows)  and  various  sizes  (columns);  h  in  the 
n/h  expression  is  the  smoothing  parameter  calculated  by  equation  (7.18). 


Normal  bimodal  distribution 


Theoretical  distributions  used  in  the  previous  two  subsections  were  unimodal.  It  can  be 
expected,  however,  that  in  reality  multimodality  is  not  an  exceptional  case.  The  results  of 
numerical  experiments  presented  in  this  and  subsequent  sections  are  to  show  how  bimodality 
of  a  theoretical  p.d.f.  is  reflected  in  its  estimated  image,  fix).  The  kernel  estimation  is 
regarded  as  the  method  which  is  very  effective  in  estimating  multimodal  p.d.f.-s  (Scott  and 
Factor  1981,  Katkovnik  1985). 

The  first  bimodal  p.d.f.  investigated  is  the  normal  bimodal  p.d.f., 

~  d’fi(x,fi^,o^  +  (\~d) 'fi(x, (”7-26) 


which  is  a  mixture  of  two  normal  density  functions, i=l,2: 


(x-ixfi 
2  a- 


(7.27) 


Parameter  d  can  be  interpreted  as  the  probability  of  taking  a  random  number  from  population 
f„(x\ixi,Oi).  For  the  purpose  of  simulation  in  all  cases  a  value  of  0.85  was  adopted  for  d,  as 
well  as  cr,  =  (Tj  =  jU;  =  1.  The  only  varying  parameter  is  H2  for  which  values  2,  4  and 
6  were  adopted  (columns  1,  2  and  3  in  Fig.  7.3,  respectively).  Varying  of  the  fHi  values 
increases  the  distance  /X2-jWi  between  the  maxima  of  the  component  p.d.f.-s  (7.27).  The 
distance  of  fx^rU’i  —  1  (column  1  in  Fig.  7.3)  is  too  small  to  reveal  these  maxima  in  the 
mixture  of  distributions;  the  theoretical  p.d.f.  (7.26)  does  not  show  any  trace  of  original 
bimodality.  Out  of  six  cases  in  column  1,  two  estimated  p.d.f.-s  are  bimodal.  It  is  interesting 
that  this  bimodality  seems  to  some  extent  to  indicate  the  underlying  original  maxima. 

If  the  distance  is  greater  than  in  the  previous  case,  bimodality  begins  to  be  visible 
in  a  subtle  way  (column  2  in  Fig.  7.3:  -  3)  although,  strictly  speaking,  the  resulting 

theoretical  p.d.f.  is  unimodal,  and  only  a  heavy  tail  occurred  on  the  right.  The  estimated 
p.d.f.-s  are  in  most  cases  bimodal  with  their  main  mode  near  the  theoretical  one 
(approximately  1).  The  second  mode  is  rather  loosely  related  to  the  second  original  mode 
(approximately  4). 

Third  column  in  Fig.  7.3  shows  the  theoretical  bimodal  p.d.f.  for  the  greatest  distance 
between  the  maxima:  jlij-Mi  =  5.  The  bimodality  of  the  theoretical  p.d.f.  is  clearly  visible. 
The  efficiency  of  the  kernel  method  is  in  this  case  very  high:  the  resemblance  of  the 
estimated  p.d.f. 
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to  the  theoretical  one  seems  to  depends  weakly  on  sample  size  n.  Almost  all  estimated  p.d.f.-s 
are  very  similar  to  their  theoretical  original  even  in  their  intensity  (especially  for  n>60).  As 
it  is  expected,  the  increase  in  n  improves  the  similarity  of  the  estimated  p.d.f.-s  to  their 
theoretical  original. 

Lognormal  bimodal  distribution 

This  distribution  is  a  mixture  of  two  lognormal  distributions  (7.24)  with  probability  d  of 
choosing  the  distribution  with  /^(x;  against  the  second  distribution  with  p.d.f. 

The  adopted  value  of  d,  0.8,  prefers  the  first  component  distribution  similarly  as  for  the 
normal  bimodal  distribution  in  the  previous  section.  Parameters  ja,,  of  the  first  component 
distribution  are  constant  and  their  values  are  1  and  0.3,  respectively,  which  corresponds  to 
the  maximum  x-value  =  cxp(^l-a^)  «  2.5.  Value  of  the  02  parameter  for  the  second 
component  distribution  is  also  constant  and  equals  0.1,  while  /nj  is  variable  and  its  value  are: 
1.5  (first  column  in  Fig.  7.4),  1.7  (second  column),  and  1.9  (third  column).  Maxima  of  that 
distribution  are  also  variable  and  their  values  are  4.4,  5.4,  and  6.6,  respectively.  Overlapping 
of  distributions  shifts  the  maxima  only  a  bit  closer  to  each  other. 

As  in  the  previous  subsection  for  normal  bimodal  p.d.f.,  in  almost  all  cases  in  Fig.  7.4 
the  estimated  p.d.f.-s  mimic  main  features  of  the  theoretical  p.d.f.  very  well.  For  as  low 
values  of  n  as  60  (and,  of  course,  for  the  higher  values)  the  maxima  are  clearly  revealed  and 
their  locations  are  in  phase  with  their  theoretical  originals.  This  concerns  especially  the  cases 
in  columns  2  and  3  where  the  theoretical  maxima  are  more  separated  than  those  in  column 
1.  Small  distance  between  original  maxima,  resulting  in  large  overlapping  of  component 
distributions,  generates  a  distinct  distortion  in  their  sample  images,  which  can  be  seen  in 
column  1  where  n  even  as  large  as  180  does  not  remove  this  discrepancy. 

In  general,  it  can  be  stated  that  the  shape  reproduction  of  the  theoretical  p.d.f.  by  a  sample 
p.d.f.  is  very  good  already  for  n  as  low  as  60.  As  expected,  the  increase  in  n  makes  the 
resemblance  of  the  theoretical  to  empirical  p.d.f.-s  larger. 


Fig.  7.4.  Lognormal  bimodal  p.d.f.-s:  theoretical  (equation  (7.28))  (heavy  line)  with  parameters  (Ati,a,),(M2.<^2)>  and 
d  shown,  and  kernel  estimated  (thin  line  with  circles)  for  samples  of  the  same  size  (rows)  and  various  sample  sizes 
n  (columns);  h  in  the  nih  expression  is  the  smoothing  parameter  calculated  by  equation  (7.18). 
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Edge  effect  case:  a  sinus  bimodal  distribution 


The  truncated  exponential  distribution  already  discussed  in  this  section  exhibits  the  strong 
edge  effect  which  substantially  distorts  the  sampled  image  of  the  true  p.d.f.  There  was  no 
reason  to  expect  that  the  distribution  of  deflection  would  smoothly  decay  at  the  limits  of  the 
variable  range.  We  could,  however,  mitigate  edge  effects  using  the  angle  identity  of 
deflection  determined  by  (6.1).  It  was  thus  important  to  know  how  weakening  of  the 
sharpness  of  the  jumps  in  a  true  p.d.f.  would  affect  estimated  p.d.f.-s.  For  this  reason  the 
following  distribution  defined  by  a  simple  p.d.f.  was  used: 


—{l+  bsinix)),  xG(c,37r-c) 
a 


0,  otherwise 


where  the  normalizing  constant  a  is  equal  to 

a  =  Stt  -  2c  +  2bcos(c) 


(7.29) 


(7.30) 


For  the  purposes  of  the  simulation  the  following  parameter  values  were  adopted:  &  =  0.8  (all 
cases),  c  =  irll  (column  1  in  Fig.  7.5),  c  =  tt/A  (column  2),  and  c  =  0  (column  3).  The 
adopted  values  of  c  determines  interesting  cases  showing  how  the  smoothness  of  the  tails  of 
a  theoretical  p.d.f.  enables  the  estimator  to  reproduce  its  original.  If  there  is  no  falling  tail 
then  the  peaks  occur  in  estimated  p.d.f.  tails  that  are  shifted  towards  the  center  part  of  the 
interval  of  realizations.  A  slightly  falling  tail  (column  2)  makes  that  shift  smaller  and  the 
greater  interval  of  the  center  part  of  the  theoretical  p.d.f.  is  better  reproduced.  In  column  3, 
locations  of  maxima  of  estimated  p.d.f.-s  are  placed  almost  at  the  same  x-values  as  the 
theoretical  maxima. 


romnarison  of  two  methods  of  estimation  of  the  smoothing  parameter  h 

Two  methods:  the  modified  (empirical)  likelihood  criterion  method  (7.22)  (or  the  cross- 
validation  method)  and  the  simplest  method  (7.18)  based  on  minimization  of  the  IMSE  of  the 
theoretical  p.d.f.  f{x)  are  compared  in  this  section.  Based  on  the  real  data  sets  of  deflections 
and  logarithms  of  energy  of  tremors,  the  comparison  is  to  show  how  both  methods  work  and 
give  some  additional  justification  for  choosing  of  the  method  (7.18). 
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The  data  sets  for  the  study  were  selected  so  that  they  represented  the  variety  of  empirical 
distribution  shapes  of  both  logarithm  of  energy  and  deflection.  The  sample  notation  used  in 
this  section  is  given  in  Table  7.1. 


Table  7.1.  Studied  data  sets 


Data  set  name 

Variable 

Origin  of  the  data 

501_CEN 

deflection 

events  from  norther  part  of  C  region  at 
level  501 

SW_E5 

deflection 

events  of  energy  >  lO^J  from  SW  region 
at  level  501 

SW_45 

deflection 

events  of  energy  >3.2"10''J  from  SW 
region  at  level  501 

416_LOGE.SW 

logarithm  of  energy 

events  from  SW  region  at  level  416 

504_LOGE.NW 

logarithm  of  energy 

events  from  NW  region  at  level  504 

510_LOGE.NEA 

logarithm  of  energy 

events  from  NE  region  at  level  510 

The  modified  likelihood  criterion  method  (7.22)  has  been  investigated  very  extensively  for 
various  kernels,  various  sample  sizes  and  various  distribution  functions  (Katkovnik  1985, 
Scott  and  Terrel  1987).  For  purposes  of  this  paper,  two  conclusions  can  be  interesting:  (i)  this 
procedure  has  turned  out  to  be  quite  effective  even  for  small  sample  sizes,  and  (ii) 
discontinuous  p.d.f.-s  are  reproduced  worse  than  the  continuous  ones.  Katkovnik  (1985) 
investigated  also  the  quality  of  the  nearest-neighbor  method.  If  unimodal  distributions  are 
taken  into  considerations,  the  method  is  comparable  to  the  modified  likelihood  criterion 
method.  However,  when  the  theoretical  p.d.f.  is  bimodal,  the  nearest-neighbor  method  quality 
is  much  worse. 

Figures  7.6  and  7.7  present  the  results  of  optimization  of  the  smoothing  parameter  h  by 
the  modified  likelihood  criterion  method  (left  column  in  both  Figures).  There  are  also  given 
nonparametric  p.d.f.-s  and  cumulative  distribution  functions  (cdf’s)  calculated  using  the 
optimized  value  of  h,  hi,  and,  as  comparison,  the  same  functions  employing  the  smoothing 
parameter  value  hi  calculated  by  equation  (7.18). 

All  lnL{h)  functions  are  multimodal,  with  the  number  of  modes  varying  from  three  to  five. 
It  is  possibly  that  shorter  optimization  steps  would  reveal  more  modes.  To  further  calculation 
the  values  maximizing  the  InLQi)  function  within  the  investigated  interval  were  chosen.  They 
are  indicated  as  hi  values  in  Figures  7.6  and  7.7.  Also  sample  data  are  shown  (right  columns 
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Fig.  7.6.  The  modified  likelihood  criterion  method  optimization  of  the  smoothing  parameter  h  for  three  real  data  sets 
of  deflections.  InlQi)  is  the  logarithm  of  the  empirical  likelihood  function  (7,21),  n  denotes  sample  size,  hjnin  is  the 
lower  boundary  value  for  optimized  h  calculated  by  equation  (7.23),  h\  and  h2  are  the  values  of  the  h  calculated  by 
equation  (7.18),  and  taken  as  the  /z-value  for  the  largest  maximum  of  lnl{h),  respectively.  Symbols/_em/7  and  F_emp 
denote  the  estimated  probability  density  function  and  the  cumulative  distribution  function,  respectively.  Scattered 
points  in  the  F_emp  plots  are  the  sample  data. 


in  Figs.  7.6  and  7.7,  the  F  emp  plot)  in  the  form  of  the  classical  empirical  distribution 
function. 

All  optimized  values  for  the  defection  p.d.f.  smoothing  parameter  (hi  in  Fig.  7.6)  are 
smaller  than  the  corresponding  simplest  method  values  hi  (for  the  sake  of  clarity,  the  values 
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Fig.  7.7.  The  modified  likelihood  criterion  method  optimization  of  the  smoothing  parameter  h  for  three  real  data  sets 
of  logarithms  of  energy.  InlQi)  is  the  logarithm  of  empirical  likelihood  fimction  (7.21),  n  denotes  sample  size,  hjnin 
is  the  lower  boundary  value  for  optimized  h  calculated  by  equation  (7.23),  hi  and  hi  are  the  values  of  the  h 
calculated  by  equation  (7.18),  and  taken  as  the  /z-value  for  the  largest  maximum  of  lnl{h),  respectively.  Symbols 
f_emp  and  F_emp  denote  the  estimated  probability  density  function  and  the  cumulative  distribution  function, 
respectively.  Vertical  heavy  lines  in  F_emp  plots  are  the  sample  data. 


of  InLQi)  for  higher  values  of  h  are  not  shown  in  Figure),  which  implies  that  using  hi  values 
will  discover  more  multimodal  structure  of  the  estimated  p.d.f.,  ^a:),  than  it  would  happen 
if  the  h\  values  were  used.  The  consequences  can  be  very  different:  from  almost  negligible, 
as  for  the  501_CEN  data  set,  to  those  introducing  great  changes,  as  it  is  for  the  SW_E45  data 
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set.  The  intermediate  case  is  shown  for  the  SW_E5  data  set. 

For  the  first  case  (501_CEN),  the  hi  value  is  approximately  a  half  of  h\.  Such  a 
difference  does  not  make  substantial  variations  to  the  estimated  p.d.f.,  although  revealing 
some  additional  rather  flat  modes. 

The  optimized  h,  hi  (equal  to  8.39),  for  the  SW_E5  data  set  is  smaller  than  that  of 
501_CEN  and  is  about  one  fourth  of  h\  (32.01).  The  resulting  p.d.f.  has  so  many  modes 
(many  of  them  are  very  subtle)  that  it  could  be  expected  that  some  of  them  do  not  reflect  the 
underlying  population  and  are  an  effect  of  sample  randomness.  The  same  effect,  in  more 
developed  version,  can  be  seen  for  the  SW_45  data  set.  It  seems  that  the  p.d.f.  using  hi  = 
3.23  reveals  much  more  noise  than  the  amount  of  population  information  that  is  hidden  by 
the  p.d.f.  using  h\  =  32.45.  All  those  relations  between  h\  and  hi  do  not  seem  to 
substantially  affect  the  cdf’s  -  all  of  them  are  little  sensitive  to  the  changes  in  the  smoothing 
parameter  values. 

For  energy  logarithms,  log{E),  Fig.  7.7,  the  general  situation  is  different  than  that  for 
deflections.  The  left  boundary,  hjnin,  of  the  range  of  optimized  values  of  h  is  in  all  cases 
greater  than  h\,  so  the  smoothness  produced  by  the  modified  likelihood  criterion  method 
(7.22)  is  greater  than  if  h\  is  used.  However,  all  hi  values  are  not  very  different  from  h\ 
values,  and  the  resulting  p.d.f.-s  are  very  similar. 

Conclusions 


1 .  Properties  of  the  nonparametric  kernel  estimation  method  of  probability  density  function 
were  investigated  in  order  to  recognize  the  possibility  of  its  application  to  objective 
reproduction  of  deflection  probability  density  functions  based  on  empirical  data.  The 
nonparametric  kernel  estimation  method,  a  novel  and  iimovative  approach  in  studies  of 
structure  of  seismic  catalogs,  is  designed  for  direct  estimation  of  an  unknown  true  probability 
function  of  a  random  variable  X  given  only  a  random  sample  x,,  i  =  l,2,...,n. 

2.  Monte  Carlo  experiments  with  variety  of  both  unimodal  and  multimodal  probability 
distributions,  and  investigations  carried  out  on  real  (empirical)  sets  of  deflections  and 
logarithms  of  energy  of  seismic  events  from  the  Wujek  mine,  made  it  possible  to  select  a 
method  of  estimation  of  the  smoothing  parameter  the  choice  of  which  is  crucial  for  the 
performance  of  the  kernel  method  for  actual  data.  The  investigation  results  show  that  when 
application  purposes  of  the  nonparametric  kernel  estimation  method  are  taken  into 
consideration,  the  most  suitable  estimation  method  of  that  parameter  is  the  method  of 
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minimization  of  the  integral  square  error.  The  selected  estimation  method  of  the  smoothing 
parameter  oversmoothes  to  some  extent  the  resulting  probability  density  function  estimator. 
For  that  reason,  in  analysis  of  real  data  not  only  distinct  modes  of  this  estimator  should  be 
traced  but  also  minor  inflexions  distorting  the  smooth  course  of  the  estimated  p.d.f.  should 
be  considered.  However,  the  chosen  method  of  smoothing  parameter  estimation  ensures  that 
the  nonparametric  kernel  estimation  method  is  independent  of  the  user  and  its  oversmoothing 
protects  against  occurrence  of  artefacts  on  the  resulting  p.d.f. 

3.  Further  on  carried  out  experiments  showed  that  for  samples  of  sizes  of  50  -  60  the 
nonparametric  kernel  estimation  method  well  reproduces  both  the  number  of  modes  of  true 
probability  distributions  and  the  locations  of  these  modes.  An  interesting  feature  of  the 
nonparametric  kernel  estimation  method  is  the  fact  that  sometimes  an  estimated  p.d.f.  enables 
better  identification  of  distribution  modality  than  the  course  (plot)  of  the  true  p.d.f. 

4.  When  the  true  p.d.f.  does  not  tend  to  zero  at  the  variability  interval  limits  of  a  random 
variable,  the  nonparametric  kernel  estimation  method  produces  an  edge  effect.  This  makes 
estimated  p.d.f. -s  smaller  than  the  true  p.d.f. -s  near  the  limits  of  the  variability  interval  and 
may  generate  false  modes  near  these  limits.  The  effect  does  not  affect  the  central  part  of  the 
estimated  p.d.f.  Because  we  expect  that  the  nonparametric  kernel  method  for  deflections  will 
be  endangered  by  the  edge  effect,  in  estimation  of  deflection  distributions  we  will  use  the 
angle  identity  of  deflection  (6.1),  estimate  p.d.f. -s  several  times  for  cyclically  transformed 
deflection  values,  and  glue  their  central  parts. 

5.  In  future,  the  nonparametric  kernel  estimation  method  can  be  used  for  estimation  of 
distributions  of  energy  or  logarithm  of  energy  of  the  seismic  events  based  on  empirical  data. 
However,  because  of  the  recently  applied  methods  of  evaluation  of  seismic  energy,  empirical 
energy  distributions  behave  as  discrete  ones.  Further  studies  are  needed  to  distinguish  the 
modes  which  are  produced  by  the  repeatability  of  some  energy  values  as  the  consequence  of 
energy  evaluation  methods  from  real  modes  which  are  the  result  of  the  multimodality  of 
energy  distribution. 


8. 


Deflections.  Identification  of  modes  of  distribution. 


The  nonparametric  kernel  estimation  of  probability  density  function  of  deflections  was 
applied  to  all  data  sets  given  in  Table  6. 1  in  Section  6. 

It  was  shown  in  Section  7,  that  when  the  actual  probability  distribution  does  not  tend  to 
zero  at  limits  of  the  variability  range  of  the  random  variable  the  nonparametric  kernel  estimation 
of  the  p.d.f.  introduces  some  boundary  errors  or  edge  effects.  To  eliminate  them  from  the 
estimated  p.d.f. -s  we  make  use  of  the  angle  identity  of  deflection  (6.1).  For  each  data  set  the 
deflections  were  shifted  four  times,  each  time  of  30° ,  so  that  the  same  data  covered  the  ranges 
(-90°,  90°),  (-60°,  120°),  (-30°,  150°)  and  (0°,  180")  respectively.  Next,  the  p.d.f.  was 

evaluated  for  each  of  the  four  ranges.  The  final  form  of  the  p.d.f.  was  obtained  by  gluing  central 
parts  of  these  four  estimates.  Finally,  the  resultant  p.d.f.  was  re-transformed  to  the  (-90°,  90°) 
range. 

The  shape  of  the  unimodal  distribution  should  be  smooth.  Every  distortion  of  the  p.d.f. 
like  high  order  discontinuous  inflexions,  folds  etc.  usually  represent  summarised  effects  of  two 
or  more  peaks.  Our  Monte-Carlo  simulations  and  the  real  data  studies  of  the  p.d.f.  kernel 
estimation  method  for  multimodal  distributions  point  out  that  the  estimates  particularly  well 
reproduce  locations  of  maxima  and  inflexions  caused  by  the  superposition  of  peaks.  The 
smoothing  parameter,  which  we  accept,  smoothes,  however,  the  final  p.d.f.  form.  Therefore, 
besides  the  distinct  peaks  of  the  estimates,  the  minor  effects  should  also  be  traced  and 
interpreted. 

As  mentioned,  it  is  expected  that  plane  structures  (linear  on  the  horizontal  plane),  like 
old  workings  edges,  faults  and  other  weak  zones  have  a  particular  significance  for  the  seismic 
event  generation.  Therefore,  we  tried  to  correlate  the  dominant  and  secondary  trends  identified 
in  the  distributions  of  deflections  with  the  directions  of  such  plane  structures  picked  firom  the 
mine  maps. 

The  results  of  the  deflection  analysis  for  the  particular  mine  parts  are  given  below. 


The  series  recorded  in  C  part,  at  the  coal  level  416 

The  estimated  p.d.f.-s  of  deflection  for  the  seismic  events  recorded  in  C  part,  at  the  level 
416  (data  sets  1,  2,  3)  are  shown  in  Figure  8.1.  In  all  distributions  the  main  mode  has  the 
maximum  located  at  the  same  place,  about  2°.  This  direction  perfectly  agrees  with  the  direction 
of  front  advance  of  aU  active  longwalls  in  this  mine  region  (see:  Appendices  1  and  2).  The  main 
mode  of  the  deflection  distribution  for  the  series  of  events  of  energy  from  lO^'j,  has  the  complex 
structure  and  the  secondary  peak  at  some  -17°  can  be  seen.  This  secondaiy  mode  can  also  be 
traced  for  the  series  with  threshold  3.2xl{f  J.  No  mine  structure  can  be  correlated  with  the 
direction  of  -17°.  The  small  peak  is  visible  in  the  distribution  for  strong  event  series  (threshold 
lO^^J),  at  some  -72°.  This  direction  can  be  correlated  with  the  strike  direction  of  Klodnicld  fault 
or  with  the  direction  of  the  edge  of  gob  in  slab  2  of  the  coal  layer  416.  Minor  distortions, 
however  presait  in  all  distributions,  are  identified  at  some  -30°  or  -35°  and  at  some  47°.  While 
the  first  direction  cannot  be  correlated  with  any  mine  structure,  the  second  one  agrees  with  the 
direction  of  the  small  fault,  of  some  2  m  throw,  located  leftside  the  seismically  active  area. 

The  series  recorded  in  SW  part,  at  the  coal  level  416 

This  is  the  first  case  of  “strange”,  that  is  not  justified  by  active  excavations,  seismicity. 
Similar  seismic  cluster  of  unknown  origin  was  observed  in  this  part,  at  the  level  507.  It  was 
speculated  that  these  clusters  could  be  side  results  of  the  tremor  generation  processes,  that  took 
place  in  this  part,  at  the  level  501.  The  estimated  p.d.f.  of  deflections  for  the  seismic  cluster 
from  the  level  416  is  given  in  Figure  8.2.  Two  distinct  modes  of  the  distribution  are  located  at 
some  -86°  and  -47°.  The  first  one  and  also  little  jumps  at  some  -68°  and  -62°  correlate  with  the 
strike  directions  of  Klodnicki  fault,  which  enclose  the  region  from  the  south  and  bends  in  this 
area  (see:  Appendices  1,  2).  The  other  one  (-47°)  cannot  be  correlated  with  any  mine  structure. 
This  mode  is,  however,  also  present  in  distributions  of  deflections  for  seismic  series  recorded  in 
the  same  part  of  the  mine,  at  the  levels  501  and  507.  The  location  of  the  small  but  distinct 
distortion,  at  about  26°,  agrees  roughly  with  the  average  of  .strke  directions  of  Arkona  (13m 
throw)  and  the  other  secondary  (4  meter  throw)  faults.  These  faults  unite  in  this  region,  and  the 
seismically  active  area  was  located  at  this  comer.  Finally  the  distortion  around  -29°  cannot  be 
linked  to  any  mine  stmcture. 
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Figure  8.2  Probability  distribution  function  of  deflection  of  events  recorded  in  SW  region, 
at  the  level  416 


The  series  recorded  in  SW  part,  at  the  coal  level  501 


The  estimated  p.d.f.-s  of  deflections  for  the  seismic  series  recorded  in  SW  part,  at  the 
level  501  (data  sets  5,  6,  7,  8)  are  shown  in  Figure  8.3.  The  identified  modes  and  their  possible 
correlations  with  mine  structures  are  as  follows: 

i)  52®.  The  maximum  of  the  main  peak  for  set  5,  splits  for  set  6,  is  hardly  visible  for  set  7  and 
practically  disappears  for  set  8.  No  structure  matches  this  direction.  It  is  likely  that  this 
maximum  value  of  the  main  mode  of  distribution  for  all  recorded  events  is  incidentally  caused 
by  summing  up  poorly  separated  mode  components; 

ii)  some  67®  to  77°.  The  distortion  in  the  distribution  for  set  5,  can  be  clearly  identified  as  a 
component  of  the  main  mode  for  set  6  and  becomes  the  main  peak  maximum  for  sets  7  and  8. 
The  direction  perfectly  agrees  with  the  strike  direction  of  Arkona  fault  of  some  20m  throw  and 
with  the  front  advance  directions  of  the  majority  of  active  longwalls  in  this  area  (see: 
Appendices  3,  4,  5).  The  increasing  significance  of  this  direction  for  the  seismic  series  from 
which  small  events  were  removed  evidences  the  particular  role  played  by  faults  in  generating 
strong  tremors; 

iii)  -22°.  The  second  distinct  mode  for  set  8.  It  is  also  visible  in  the  distribution  for  set  7  but  its 
presence  in  the  distribution  for  sets  6  and  5  can  only  be  guessed.  This  is  the  strike  direction  of 
all  longwalls  worked  in  the  area.  The  open  question  remains  why  this  direction  is  important  only 
for  the  strong  event  series; 

iv)  some  -5®  or  -8®.  Appear  as  distortions  in  sets  5  and  6  curves  and  as  a  distinct  maximum  in  set 
7  curve.  It  is  not  present  in  the  distribution  for  set  8.  No  simple  geometrical  explanation  can  be 
given  for  this  mode; 

v)  some  22°  to  27®.  The  distortions  of  the  p.d.f.  estimates  at  these  angles  are  observed  for  all 
samples.  The  direction  agrees  with  the  strike  direction  of  the  second  fault  (4.5m  throw)  in  the 
area; 

vi)  some  minor  effects  appear  here  and  there  at  -70®  or  -80®  and  at  -40°  or  -45®. 

One  can  note  that  the  dominant  direction  of  the  series  rotated  of  as  much  as  20°  after 
small  events  were  removed  from  the  sample.  Besides  changing  location  of  maximum,  the  main 
peak  of  the  p.d.f.-s  of  deflections  became  narrower  while  increasing  the  energy  threshold. 


Moreover,  the  large  tremor  series  have  also  the  distinct  secondary  trend  which  is  not  visible  in 
the  distribution  for  the  series  comprising  all  recorded  events. 

The  series  recorded  in  the  southern  part  of  C  region,  at  the  coal  level  501 

The  estimated  p.d.f.-s  of  deflections  for  to  the  data  sets  9,  10,  11,  12  are  shown 
respectively  in  Figure  8.4.  One  may  see  the  main  mode  of  the  distribution  for  the  series  of  aU 
recorded  events  (set  9)  located  in  the  range  from  80°  to  -60°.  the  mode  narrows  and  shifts  its 
maximum  from  -82°  to  -1-81°  after  the  threshold  of  lO^^J  is  applied  and  uncovers  a  complex 
subtle  structure  in  the  range  from  -75°  to  -50°.  These  two  new  components  again  unite  after 
increasing  the  threshold  of  the  next  half  of  order  (set  11).  Finally,  for  the  series  comprising  the 
strongest  events  (set  12)  the  main  mode  of  the  distribution  of  deflections  is  narrow  and  regular, 
with  maximum  at  90° .  90°  is  the  strike  direction  of  all  active  longwaUs  in  C  part,  at  this  coal 
level  (see:  Appendices  3, 4,  5). 

Other  distortions  of  the  set  9  curve  are  visible  at  some  15°  and  from  0°  to  -20°  .  The  are 
not  displayed  in  the  p.d.f.  for  set  10  except  a  slight  trace  at  some  15°  but  they  reappear  in  the 
distribution  for  set  11.  Finally  the  distinct  and  complex  secondary  mode  is  formed  in  the  range 
from  -20°  to  20°  for  the  series  of  strong  events  (set  12).  This  mode  can  be  correlated  with  the 
front  movement  directions  of  all  longwaUs  and  the  direction  of  old  gob  edges,  in  the  area. 

The  series  recorded  in  the  northern  part  of  C  region,  at  the  coal  level  501 

Set  13  is  the  other  "strange"  case  among  induced  seismicity  clusters  from  Wujek  mine. 
As  mentioned  previously,  no  miming  works  were  carried  on  in  this  part  of  the  mine  during  the 
period  under  study.  The  estimated  p.d.f.  of  deflections  for  this  series  of  seismic  events  is  shown 
in  Figure  8.5.  The  main  mode  is.wide  with  maximum  at  some  -70°  .  This  mode  could  be 
eventually  correlated  with  the  closest  edge  of  old  gobs  in  slab  3  of  this  coal  layer  (see:  Appendix 
3).  Minor  distortions  of  the  p.d.f.  curve  are  seen  at  some  -55°,  -f-46°  and  -f  15°  respectively.  No 
mine  structures  can  be  related  to  these  directions. 
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Figure  8.5  Probability  distribution  function  of  deflection  of  events  recorded  in  the  northern 
part  of  C  region,  at  the  level  501 
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The  series  recorded  in  C  part,  at  the  coal  level  504 

The  estimated  p.d.f.-s  of  deflections  for  the  seismic  series  recorded  in  C  part,  at  the  level 
504  (data  sets  14,  15,  16)  are  shown  in  Figure  8.6.  Besides  the  main  mode  with  maximum  at 
about  l°or  -5°,  present  in  all  curves,  one  can  observe  development  of  the  second,  wide  mode  at 
some  55°or  60°  with  the  increase  the  energy  threshold  value  applied  to  the  data.  For  the  set 
comprising  the  strongest  event  this  second  mode  is  nearly  as  big  as  the  main  one.  The  main 
mode  direction  matches  the  direction  of  advance  of  all  active  longwalls  in  this  part  of  the  mine 
(see;  Appendix  6).  The  second  mode  direction  agrees  with  the  strike  direction  of  the  fault  of 
some  2  m  throw  which  crosses  the  area  of  excavations.  The  figure  gives  the  clear  evidence  of 
fault  influence  on  generating  of  strong  events. 

There  are  other  distortions  of  the  deflection  p.d.f.-s  at  some  -18°  or  -20°  ,  35°  or  40°  and 
the  most  distinct  one  between  -65°  and  -75°.  The  last  direction  corresponds  to  the  strike  direction 
of  the  Klodnicki  fault,  that  encloses  the  area  from  the  south.  For  the  other  directions  no  clear 
link  to  mine  structures  can  be  found. 

The  series  recorded  in  NW  part,  at  the  coal  level  507 

The  estimated  p.d.f.-s  of  deflection  for  the  data  sets  17,  18,  19  are  shown  in  Figure  8.7. 
The  main  mode,  with  maximum  located  at  some  -80°,  becomes  wider  and  displayes  a  complex 
multicomponent  structure  after  the  prescribed  energy  thresholds  are  applied  to  the  data.  One  can 
see  the  second  peak  in  the  distribution  of  deflections  for  the  strongest  event  series  at  75°.  The 
-80°  direction  is  the  direction  advance  of  the  active  longwalls  in  this  area  (see:  Appendix  7).  The 
75°  direction  matches  the  strike  direction  of  the  nearby  Arkona  fault  (15-25  m  throw  in  this 
area).  Again  the  influence  of  the  fault  becomes  visible  when  only  strong  events  are  taken  into 
account.  The  slight  jump,  located  at  4°  by  the  global  minimum  of  distribution  for  the  strong 
event  series,  correlates  with  the  strike  direction  of  the  small  fault  wing  (1.8  meter  throw)  that 
crosses  the  area  of  excavation.  The  distortion  noticed  in  the  same  curve  at  some  53°  caimot  be 
linked  to  mine  structures. 
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The  series  recorded  in  SW  part,  at  the  coal  level  507 

Here  we  have  the  next  example  of  the  seismic  cluster  that  occurred  in  the  area  where  no 
mining  works  were  carried  on.  As  mentioned,  the  similar  effect  was  observed  in  the  same  part 
of  the  mine,  at  the  level  416.  The  estimated  p.d.f.  of  deflections  for  this  cluster  (set  20)  is 
presented  in  Figure  8.8.  The  main  mode  of  the  distribution  has  comlex  shape  and  is  located 
between  -48°  and  -38°  .  There  is  no  structural  explanation  for  this  direction.  However,  the  same 
direction  modes  have  been  also  found  for  events  recorded  in  the  same  region,  at  the  level  416 
(the  main  mode)  and  at  the  level  501  (the  minor  effect).  The  other  identified  distortions  are: 

i)  The  flat  well  separated  peak  at  some  14°.  The  direction  agrees  with  the  strike  direction  of  the 
small  fault  (3-4m  throw)  located  rightside  of  the  area  (see:  Appendix  7); 

ii)  Some  74°.  This  direction  matches  the  strike  direction  of  Arkona  fault,  located  leftside  of  the 
area.  It  also  agrees  with  the  dominant  trend  of  the  seismic  series  recorded  in  this  area,  at  the 
level  501; 

iii)  possibly  about  -73°.  No  mine  structure  can  be  correlated  with  this  eventual  mode. 

The  series  recorded  in  NW  part,  at  the  coal  level  510 

The  estimated  p.d.f. -s  of  deflections  for  the  seismic  series  recorded  in  NW  part,  at  the 
level  510  (data  sets  21,  22,  23)  are  presented  in  Figure  8.9.  The  case  turns  down  eventual 
suspitions  that  the  rising  complexity  of  the  distributions  of  deflections  with  the  increase  of  the 
value  of  energy  threshold  value  could  be  caused  by  the  decrease  of  sample  sizes.  Unlike  in  the 
majority  of  cases,  here  the  shape  of  the  estimated  p.d.f.  of  the  seismic  series  composed  of  strong 
events  (set  23)  is  less  complex  than  that  for  the  series  comprising  all  recorded  events  (set  21). 
The  main  mode  of  the  distribution  for  the  series  of  all  events  has  the  maximum  located 
somewhere  between  -85°  and  -65°.  For  the  strong  event  series  the  maximum  of  the  main  mode 
of  distribution  is  identified  at  some  -75°  or  -68°.  This  direction  agrees  with  the  direction  of  front 
advance  of  the  active  longwalls  (see:  Appendices  8,  9).  Neither  the  peak  at  about  -25°  nor  the 
peak  at  some  5°,  that  are  identified  for  the  all  event  series  distribution  and  that  gradually  disapear 
after  the  prescribed  thresholds  are  appUed,  can  be  linked  to  the  mine  structures.  Two  other 
distortions  appear  in  the  distribution  for  the  series  of  events  >  lO^^J,  at  some  52°  and  possibly  at 
about  20°.  The  first  one  cannot  be  correlated  with  any  mine  structure,  but  it  also  appeared,  as 
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Figure  8.8  Probability  distribution  function  of  deflection  for  events  recorded  in  SW  region, 
I  at  the  level  507 
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Figure  8,9  Probability  distribution  functions  of  deflection  for  events  recorded  in  NW  region, 
at  the  level  510 
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the  slight  distortion,  for  strong  event  sequence  recorded  in  the  same  part  of  the  mine,  at  the  level 
507.  The  second  one  could  be  eventually  related  to  the  small  fault  adjoining  the  area  of 
excavation  from  the  right  hand  side. 

The  series  recorded  in  the  western  part  of  ME  region,  at  the  coal  level  510 

The  last  case  of  “strange”  cluster  of  seismicity.  In  this  case  the  number  of  events  and  the 
energy  distribution  of  events  allow  also  to  study  the  series  from  which  small  events  were 
removed.  The  estimated  p.d.f.-s  of  deflections  for  these  series  (data  sets  24,  25,  26)  are  shown 
in  Figure  8.10.  The  shapes  of  the  p.d.f.-s,  particularily  the  locations  of  the  global  miniminri 
significantly  differ  between  the  distributions  for  the  all  event  series  and  for  the  strong  event 
series.  The  main  mode  maximum  follows,  more  or  less,  the  same  direction  about  -80°or  -85°  for 
all  curves.  This  direction  can  be  correlated  with  the  direction  of  the  closest  gob  edge  located 
south  from  the  seismically  active  area  (see:  Appendices  8,  9).  The  other  peaks,  identified  at  -27° 
,  -47°  and  -85°  in  the  distributions  for  sets  24  and  25;  at  -9°  in  the  distribution  for  set  24,  and  at 
25°,  61°  and  85°  in  the  distribution  for  set  26  cannot  be  correlated  with  mine  structures. 

The  series  recorded  in  SE  part,  at  the  coal  level  510 

The  last  considered  case  deals  with  the  seismic  cluster  located  in  SE  part  of  the  level 
510.  This  is  the  only  analysed  cluster  that  occurred  in  SE  part  of  the  mine.  The  estimated  p.d.f.- 
s  of  deflections  are  shown  in  Figure  8.11.  The  main  mode  is  located  similarily  in  all 
distributions.  Its  maximum,  at  some  45°  to  50°,  perfectly  agrees  with  the  advance  direction  of 
the  longwalls  lA,  DA  and  of  the  longwall  DLA  in  its  first  part,  in  slab  2,  with  the  old  works  edge 
in  slab  1  and  with  the  strike  direction  of  Wojciech  fault  (see:  Appendices  8,  9).  The  longwall 
IDA  was  bend  in  its  second  part.  One  can  see  in  Figure  8.11  a  slow  development  of  the 
secondary  mode  at  some  -5°,  5°  after  the  thersholds  are  applied  to  the  data.  This  secondary  mode 
direction  matches  the  advance  direction  of  the  second  part  of  longwall  IDA  and  the  advance 
directions  of  all  langwalls  in  slab  1  of  the  level.  The  small  distortions  appearing  here  and  there 
in  the  distributions  in  the  range  between  -35°and  -55°  could  be  roughly  correlated  with  the  strike 
direction  of  the  lonwalls  in  slab  2.  The  small  distinct  peak  at  -  73°  in  the  estimated  p.d.f.  of 


deflections  of  events  >  3.2xlO'*J  agrees  in  direction  with  the  direction  of  the  Klodnicki  fault 
strike  in  this  area. 

Discussion 

Within  the  variety  of  considered  cases  it  is  not  easy  to  find  common  rules  governing  the 
processes  of  the  seismic  events  generation  at  the  identified  dominant  directions.  Seismic  series 
from  each  part  of  the  mine  have  their  own  dominant  directions  and  the  own  ways  in  which  the 
dominant  directions  change  when  different  thresholds  are  applied  to  the  data.  Certainly  at  least 
one  distinct  mode  of  the  distributions  of  deflections  appeares  in  all  cases.  The  main  mode  usually 
becomes  narrower  with  the  increase  of  the  value  of  threshold.  The  maximum  of  the  main  mode 
is  usually  located  at,  more  or  less,  the  same  place  for  all  seismic  series  associated  with  the  same 
area  of  the  mine,  regardless  the  threshold  applied.  There  are,  however,  two  cases  where  the 
locations  of  the  global  maximum  of  deflection  distributions  differ  strongly  between  the 
distributions  for  all  events  and  for  strong  events.  Both  cases  were  noted  for  the  seismic  clusters 
associated  with  the  level  501,  the  most  seismically  active  and  the  most  complex  level  from  the 
mining  point  of  view. 

The  number  of  visible  modes  in  the  distributions  of  deflections  for  the  series  of  strong 
events  is,  in  general,  greater  than  that  in  the  distributions  for  all  event  series.  One  could  suspect 
that  the  new  modes  foimd  in  the  distributions  for  the  strong  event  series  are  artefacts,  caused  by 
the  drop  of  the  sample  sizes  after  applying  the  thresholds.  There  are,  however,  cases  for  which 
the  proportion  of  the  numbers  of  modes  for  the  all  event  series  and  for  the  strong  event  ones  is 
reversed  (e.g.,  the  seismic  series  from  NW  part,  at  the  level  510).  Thus  we  think  that  this 
general  feature  is  indirectly  resulted  by  the  decreasing  complexity  of  the  series  with  the  rise  of 
the  threshold.  In  the  series  comprising  all  events  recorded  in  the  area  small  events  outnumber 
large  ones.  The  dominant  trends  of  the  small  event  deflections  change  in  time  because  of 
variations  of  the  mining  process.  Thus  in  such  series  from  long  time  periods  these  varying  trends 
construct  a  sort  of  random  noise  which  is  added  to  the  constant  components  of  the  distribution  of 
deflections.  The  structure  of  the  distribution  is  usually  so  dominated  by  this  random  noise  that 
the  secondary  modes  are  hidden.  The  tendecy  to  narrow  the  main  peak  after  small  events  are 
removed  from  the  sample  and  the  presence,  in  some  cases,  of  uninterpretable  directions  for  the 
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series  of  all  events,  which  rotate,  split  and  match  characteristic  directions  of  mine  geometry  after 
the  threshold  is  applied,  support  this  reasoning. 

In  all  but  one  considered  cases  of  seismic  activity  occurring  in  the  areas  of  active 
excavations  the  locations  of  the  main  modes  of  distributions  of  deflections  agree  with  the 
advance  directions  of  longwalls.  The  agreement  is  usually  better  for  the  distributions  for  the 
strong  event  series  than  for  those  for  the  aU  event  sets.  There  is  one  case  for  which  the  longwall 
faces  in  two  slabs  of  the  same  area  were  advancing  in  different  directions  (SE  part  at  the  coal 
level  510).  In  tMs  case  the  distribution  of  deflection  for  the  strong  event  series  displayes  two 
dominant  trends  that  respectively  match  the  advance  directions. 

For  the  seismically  active  areas  of  mining  works  at  the  level  501  the  distributions  of 
deflections  display  dominant  trends  corresponding  with  the  strike  directions  of  the  longwalls. 
The  feature  is  not  found  for  the  seismic  series  from  other  levels. 

The  distinct  secondary  directions  of  the  seismic  series  recorded  in  the  areas  of  mining 
works  usually  correspond  with  the  strike  directions  of  the  faults  located  close  to  the  areas.  This 
correspondence  is  better  for  the  strong  event  series  than  for  the  small  event  ones  evidencing  the 
importance  of  faults  in  the  strong  event  generation  processes.  Sometimes  the  secondary  mode 
directions  correlate  with  the  positions  of  edges  of  old  works.  There  are,  however,  modes  of  the 
distributions  of  deflections,  which  carmot  be  linked  neither  with  any  tectonic  nor  with  any  mine 
structures. 

For  the  seismic  clusters  that  appeared  apart  from  the  mine  work  areas  the  dominant 
directions,  found  in  the  distributions  of  deflections,  usually  correlate  with  the  strike  directions  of 
the  faults,  if  the  faults  are  present  in  the  area,  and  sometimes  with  the  directions  of  closest  edges 
of  old  works.  Again,  however,  one  can  see  modes  of  the  distributions,  which  caimot  be 
explained  by  the  locations  of  mine  structures.  Such  the  interesting  cases  concern  the  seismic 
clusters  associated  with  SW  parts  of  the  levels  416  and  507.  For  both  cases,  the  main  modes  are 
located  at  the  same  angles.  Moreover,  the  same  angle  is  assigned  by  the  secondary  mode  of 
distribution  of  deflections  for  the  series  recorded  in  this  area  at  the  level  501,  however,  no  mine 
structure  justify  importance  of  this  angle. 
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Conclusions 


1.  Directional  structure  of  the  seismic  series  is  assessed  by  tracing  the  dominant  and 
secondary  directions  in  the  probability  density  functions  of  deflections,  estimated  by  the 
nonparametric  kernel  estimation  method. 

2.  The  results  of  the  analysis  of  variuos  seismic  clusters  from  Wujek  coal  mine  show  that 
the  epicentre  distribution  of  seismic  events  there  is  not  random.  The  dominant  directions  of  the 
distribution  are  significantly  different  for  the  studied  sets  of  events  that  have  occurred  in 
different  areas  of  the  mine. 

3.  Furthermore,  the  dominant  trends  in  the  same  set  of  events  become  different  when 
various  energy  thresholds  are  considered.  Several  secondary  trends  are  also  observed  in  the 
distributions.  This  confirms  that  the  generation  of  events  is  controlled  by  various  generating 
processes  and  that  mining-induced  seismicity  is  a  multimodal  phenomenon. 

4.  It  is  possible,  by  careful  inspection  of  the  distributions  of  deflections,  to  identify 
generating  processes  of  seismicity  responsible  for  particular  modes  of  the  distribution. 

5.  In  some  cases,  after  complementing  the  results  of  the  analysis  of  deflections  with 
information  from  other  sources,  for  instance  from  mine  maps,  it  is  also  possible  to  conclude 
about  factors  that  control!  these  processes. 

6.  The  dominant  trends  of  the  seismic  series  from  Wujek  coal  mine  usually  correlate  with 
the  directions  of  front  advance  of  working  faces. 

7.  Some  of  the  secondary  trends  can  be  associated  with  known  geological,  tectonic  and 
mining  features,  whereas  the  others  are  not  directly  connected  with  the  known  mine  structures. 

8.  The  study  evidences,  in  particular,  that  faults  have  a  significant  role  for  the  strong 
seismic  event  generation.  Their  presence  in  or  in  vicinity  of  the  seismically  active  areas  changes 
the  distributions  of  deflections. 

9.  The  analysis  of  deflection  forms  a  basis  for  separating,  from  the  seismic  series,  events 
generated  by  indovidual  processes.  These  will  be  the  events  connected  with  particular  modes  of 
the  distribution  of  deflections. 

10.  The  dominant  trends  considered  here  are  based  on  bilateral  movements  of  seismic 
events.  Their  possible  migration  (for  which  unilateral  movements  must  be  considered),  a  well- 
studied  phenomen  from  natural  earthquakes  (e.g.  Gedney  et  al.  1980,  Meyer  et  al.  1985,  King 


—  Ill  — 


and  Ma  1988,  Yoshida  1988,  Rydelek  et  al.  1990,  Loo  et  al.  1992),  is  another  highly  interesting 
subject  awaiting  future  studies.  Correlations  of  the  dominant  trends  with  source  mechanism  and 
source  parameters  are  also  needed  for  better  understanding  of  the  multimodal  character  of 
seismic  events  induced  by  mining. 
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9.  Separating  of  events  generated  by  tbe  single  process 


The  analysis  of  deflections  presented  in  Sections  6,  7,  8  can  be  used  to  separate,  from  the 
seismic  series,  events  generated  by  the  processes  related  to  particular  modes  of  the  distribution 
of  deflections.  Here  one  example  of  such  the  separation  is  presented. 

The  seismic  series  recorded  in  C  part,  at  the  level  504  was  considered.  The  analysis  of 
deflection  of  the  series  showed  apparent  development  of  the  second  dominant  trend  in  the 
distribution  of  deflections  after  increasing  energy  threshold  was  applied  to  the  data  (see:  Figure 
8.6).  This  trend  was  linked  to  the  small  fault  that  crossed  the  seismically  active  area.  The  other 
dominant  trend  matched  the  direction  of  front  advance  of  all  mine  longwalls  that  were  worked  in 
the  studied  area. 

The  separation  procedure  starts  form  assigning  the  interval  of  deflections  which 
contributed  to  the  trend  of  interest.  In  this  case,  on  the  basis  of  the  estimated  p.d.f.  of 
deflections  for  the  series  of  events  of  energy  >  lO^'j,  the  interval  [45®,  75°]  was  accepted.  A 
value  of  deflection  is  evaluated  from  locations  of  epicenters  of  two  consecutive  events.  Linking 
it  with  one  of  these  events  is  artificial  and  can  only  be  used  to  study  serial  properties  of 
deflections.  In  the  separation  of  events  presented  here  we  mark  both  events  of  each  pair  in  the 
studied  series  of  events,  which  have  deflection  values  in  the  specified  interval.  According  to  our 
reasoning  given  in  previous  section  the  marked  events  are  likely  to  be  outcomes  of  the  same 
identified  generating  process.  It  was  interesting,  however  not  unexpectable,  that  the  marked 
events  form,  in  some  cases  sequences  of  more  than  two  events. 

The  spatial  distribution  of  all  seismic  events  recorded  in  the  studied  area  of  the  mine  is 
presented  in  Figure  9.1.  Positions  of  the  separated  events  are  shown  by  crosses.  The  strike  of  the 
fault,  which  could  be  the  reason  of  the  considered  dominant  trend  of  the  series,  is  approximately 
sketched  with  the  broken  line.  The  full  line  shows  boundaries  of  the  excavations  worked  during 
the  period  selected  for  our  study.  The  longwalls  were  moving  towards  the  north.  The  separated 
events  locate  on  lines  more  or  less  parallel  to  the  strike  of  the  fault.  It  is  interesting  that  they  are 
not  placed  at  the  fault  but  are  distanced,  in  average,  more  than  100m  in  the  northernmost  part  of 
the  area.  Similar  behavior  of  distributions  of  sources  were  noticed  for  other  series  of  events 
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Figure  9.1  Distribution  of  seismic  sources  in  C  part,  at  the  level  504.  Crosses  mark  event! 
separated  by  the  analysis  of  deflections.  Boundaries  of  active  excavations  are 
sketched  with  full  line  and  the  arrow  marks  the  direction  of  front  movements. 
Broken  line  marks  approximately  the  strike  of  the  fault. 


determined  as  related  to  faults.  In  some  cases  they  were  located  on  more  than  one  lines  parallel 
to  the  strikes  of  the  faults  on  the  same  fault  wings. 
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10.  Conclusions 


1 .  Selected  statistical  methods  are  proposed  to  investigate  complexity  of  the  series  of 
seismic  events  with  the  aim  to  identify  different  event  generating  processes.  The  pseudo¬ 
multivariate  approach,  that  is  the  repeated  analysis  of  one  parameter  distributions  for  the 
subsets  of  the  seismic  catalog,  extracted  according  to  conditions  put  on  the  other  parameters 
of  events,  is  used. 

2.  The  seismic  catalog  from  Wujek  coal  mine.  Upper  Silesian  Coal  Basin,  Poland 
satisfied  initial  requirements  concerning  size,  quality,  variability  and  completeness  of  the 
database.  Altogether  14297  seismic  events  were  analyzed.  Inspection  of  the  catalog 
evidences  completeness  of  the  database  and  allows  to  recognize  qualitatively  its  features, 
important  for  its  future  division  into  more  homogeneous  subseries. 

3.  The  descriptive  statistics  studies  show  that  the  induced  seismicity  is  a  complex, 
multimodal,  both  time  and  space  dependent  phenomenon.  The  event  rate  and  the  energy 
distribution  of  events  are  time  varying  due  to  time  changes  of  mining. 

4.  Moreover,  the  descriptive  statistics  studies  show  that  the  multimodality  of  the 
Seismic  series  is  reflected  in  the  distribution  of  event  energy.  At  least  two  classes  of 
generating  processes  have  been  identified.  The  one  comprises  time-varying,  controlled  by 
mining  works  processes,  responsible  for  occurrence  of,  in  general,  small  events.  The  other 
one  comprises  processes  weakly  dependent  on  time,  giving  rise  to  strong  events.  The 
processes  from  the  first  class  are  predominantly  dependent  upon  the  coal  layer  and  stope 
characteristics.  The  second  class  processes  depend  more  upon  joint  influence  of  structures  in 
the  environment  of  excavations  than  upon  the  excavations  themselves. 

5.  The  maximum  energy  value  statistics  method  complemented  with  testing 
procedures  assessing  differences  in  the  empirical  distributions  of  maxima  turns  out  to  be  a 
valuable  tool  in  the  detailed  studies  of  the  multimodal  structure  of  seismic  series.  The 
results  of  the  maximum  value  analysis  prove  that  the  energy  distribution  of  seismic  events  is 
multimodal  in  both  the  mine-wide  scale  as  well  as  for  particular  clusters  of  seismicity.  The 
method  enables  to  single  out  individual  modes  of  the  event  energy  distributions  and  to 
estimate  their  limiting  values.  However,  the  extreme  value  statistics  method  is  an  indirect 
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way  to  identify  processes  of  seismic  event  generation  and  it  is  not  possible,  on  the  basis  of 
this  method,  to  separate,  from  the  seismic  series,  events  generated  by  the  specific  identified 
processes. 

6.  A  new  method  is  proposed  to  study  nonrandom  character  of  the  spatial 
distribution  of  seismic  series.  The  method  is  based  on  deflection  of  straight  line  connecting 
epicenters  of  eveiy  two  consecutive  events,  measured  from  the  N-S  direction. 

7.  The  efficiency  an  the  conditions  of  applicability  of  the  nonparametric  kernel 
estimation  method  of  probability  density  function  for  reproducing  of  deflection  probability 
functions,  based  on  empirical  data,  have  been  investigated.  The  nonparametric  kernel 
estimation  method,  a  novel  and  innovative  approach  in  studies  of  structure  of  seismic 
catalog,  is  an  objective,  user-independent  technique  which  can  be  used  for  deflection  series 
of  sizes  from  50-60.  The  method  well  reproduces  both  the  number  of  modes  of  true 
probability  distributions  and  the  locations  of  these  modes. 

8.  In  future  the  nonparametric  kernel  estimation  can  be  also  used  for  estimation  of 
energy  or  logarithm  of  energy  of  the  seismic  events,  based  on  empirical  data.  Further 
studies  are  needed  to  account  for  distortions  of  the  estimated  probability  density  function  of 
energy  introduced  by  the  repeatability  of  some  energy  values. 

9.  The  application  of  the  analysis  of  deflection  with  the  use  of  the  nonparametric 
kernel  estimation  of  the  probability  distribution  function  to  the  data  from  Wujek  coal  mine 
shows  that  the  epicenter  distribution  of  seismic  events  there  is  not  random.  The  dominant 
directions  of  the  distribution  are  significantly  different  for  the  studied  subsets  of  events  that 
occurred  in  different  areas  of  the  mine.  Furthermore,  the  dominant  trend,  in  the  same 
subset  of  events,  become  different  when  various  energy  thresholds  are  considered.  Several 
secondary  trend  are  also  observed  in  the  distributions. 

10.  The  main  trend  and  the  secondary  trends  as  well,  found  in  the  estimated 
probability  density  functions  of  deflections  are  supposed  to  identify  different  processes  of 
seismic  event  generation. 

11.  In  some  cases,  after  complementing  the  results  of  the  analysis  of  deflection  with 
information  from  other  sources,  it  is  possible  to  conclude  about  factors  that  control  these 
identified  processes.  Some  of  the  dominant  and  secondary  trends  in  the  distributions  of 


deflections  can  be  associated  with  known  geologic,  tectonic  and  mining  feature,  whereas  the 
others  are  not  directly  connected  with  the  known  mine  structures. 

12.  The  analysis  of  deflection  is  a  direct  method  to  study  the  seisnuc  catalog  and 
allows  to  separate,  from  the  seismic  series,  events  originated  by  the  generation  processes 
identified  in  the  distribution  of  deflections. 

13.  The  statistical  analysis  of  seismic  series  reveals  information  about  properties  of 
seismic  event  parameters  that  is  contained  in  the  data.  Its  integration,  whenever  possible, 
with  other  kinds  of  seismic  event  analysis,  e.g.  source  mechanism  and  source  parameters 
studies,  would  increase  certainty  of  identification  of  event  generating  processes  and  of 
separation  of  events  generated  by  individual  processes. 
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